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1. Introduction and History 

The 'pair-contact process with diffusion' (PCPD) is a classical reaction-diffusion 
process describing stochastically interacting particles of a single species which react 
spontaneously whenever two of them come into contact. In its simplest version the 
PCPD involves two competing reactions, namely 



In addition individual particles are allowed to diffuse. Moreover, depending on the 
specific variant under consideration, an exclusion principle or a similar mechanism may 
impose a restriction such that the particle density cannot diverge. 

During the past four years the PCPD has attracted a lot of attention. The interest 
in this model is first motivated by the fact that the PCPD exhibits a non-equilibrium 
phase transition caused by the competing character of fission and annihilation. In the 
so-called active phase, the fission process dominates and the system evolves towards a 
fluctuating steady-state characterized by a certain stationary particle density p s > 0. 
Contrarily, in the absorbing phase the annihilation process dominates and the density 
of particles decreases continuously until the system reaches a so-called absorbing state 
(for example the empty lattice) from where it cannot escape. Since absorbing states 
can only be reached but not be left, the process is generically out of equilibrium. The 
active and the absorbing phase are separated by a continuous phase transition with a 
non-trivial critical behaviour. 

As we shall see, despite the simplicity of the PCPD, the properties of its phase 
transition(s) are quite involved and rather surprisingly, at the time of writing no clear 
consensus on its critical behaviour has been achieved. This makes the PCPD an 
ideal laboratory to explore the most advanced techniques from either field-theory or 
simulation. On the other hand, due to the intensive recent research, the PCPD offers 
the possibility to introduce all these techniques as applied to the same model and thereby 
to critically appreciate their respective strengths and weaknesses. 

The classification of non-equilibrium phase transitions from an active phase into 
an inactive phase with absorbing states is a challenging problem of contemporary 
statistical physics [1,2]. Generally it is believed that such transitions belong to a 
finite number of universality classes with steady-state properties being characterized 
by a quartet f3', v±, u\\) of critical exponents and certain universal scaling functions. § 

§ Only if a certain duality symmetry holds, the a priori different [3] steady-state exponent (3 of the 



fission: 



annihilation: 



2A -> 3A , 
2A -> 0. 



Pair- contact process with diffusion 



4 



So far only few universality classes are firmly established, the most important ones 
being directed percolation (DP) [4-7], the parity-conserving (PC) class of branching- 
annihilating random walks [8,9], ^-symmetric transitions in the voter model [10,11], the 
general epidemic process [12,13], and absorbing phase transitions coupled to a conserved 
field [14-16]. The current interest in the PCPD is motivated by the perspective that 
it may represent yet another independent universality class of non-equilibrium phase 
transitions which had not been studied before. The asymptotic critical behaviour, 
however, seems to be masked by strong corrections to scaling so that it is not yet 
entirely clear how the transition in the PCPD should be characterized and classified. 
The purpose of this review is to summarize the present state of knowledge and to point 
out open questions. 

The first investigation of this class of models considered the non-diffusive pair- 
contact process (PCP) [17,18] with at most one particle per site. In this case individual 
particles are completely immobile so that activity can only spread in space by offspring 
production at empty nearest-neighbour sites. On a one-dimensional lattice, the PCP 
may be implemented by fission fyAA/AAfy — > AAA and annihilation AA — > 00. A 
typical space-time plot of a critical PCP is shown in the left panel of figure 1. Since 
solitary particles do not diffuse, any configuration of spatially separated particles is 
frozen. Therefore, on an infinite lattice the model has infinitely many absorbing states. 
More specifically, on a one-dimensional lattice with L sites and either periodic or free 
boundary conditions the model has 2 L states of which 



Numerical simulations in both one dimension (ID) [17,20] and in 2D [21] have shown 
that the steady-state transition in the non-diffusive PCP belongs to DP. Although a 
rigorous proof is not yet available, recent numerical high-precision studies of universal 
scaling functions in various dimensions support this claim [22]. The result is remarkable 
since the PCP does not obey the assumptions of the so-called DP-conjecture [23,24], 
which states that under certain generic conditions any non-equilibrium phase transition 
into a single absorbing state should belong to DP. Spreading from a single seed, however, 
may show non-universal properties depending on the density and correlations of the 
surrounding particles [3,25,26]. To explain the observed critical behaviour in the PCP 
several field-theoretic approaches have been suggested [27,28]. 

In the PCPD, where in addition to the binary PCP reactions individual particles are 
also allowed to diffuse, the situation is fundamentally different. Unlike the PCP without 
particle diffusion, the PCPD has only two absorbing states, N stat = 2. These are the 
empty lattice and a homogeneous state with a single diffusing particle. Therefore, and 

order parameter and /?' of the survival probability will agree. 




(1) 



are absorbing [19], where g + = (l + v / 5)/2 ~ 1.618 is the golden mean (see appendix A). 
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PCPD 



Figure 1. Onc-dimensional pair-contact process starting with a fully occupied lattice 
at criticality. Left: In the pair-contact process without diffusion (PCP) solitary 
particles cannot diffuse, leading to frozen patterns of separated vertical lines. Right: 
In the PCPD, where individual particles are allowed to diffuse, offspring production 
can be restarted after long times when two diffusing particles meet, leading to a very 
different visual appearance of the process. 



in contrast to spreading processes such as DP, where diffusion does not alter the critical 
behaviour [29] , introducing diffusion in the pair-contact process should be regarded as a 
singular perturbation which changes the structure of the absorbing phase. In fact, even 
the visual appearance changes drastically, as shown in the right panel of figure 1. 

The pair-contact process with diffusion was already suggested in 1982 by 
Grassberger [30], who expected a critical behaviour "distinctly different" from DP. || 
15 years later the problem was rediscovered, under the name of "annihilation-fission 
process", by Howard and Tauber [31], who proposed a bosonic field-theory for the one- 
dimensional PCPD. Because of the bosonic nature of their model, the particle density 
is unrestricted and thus diverges in the active phase. Furthermore, since their field- 
theory turned out to be non-renormalizable, no quantitative statements on the phase 
transition between the active and the inactive phase could be made. This triggered 
one of us (M.H.), together with J.F.F. Mendes, to attempt a non-perturbative study 
of the ID PCPD, formulated on a discrete lattice and with at most one particle per 
site. Starting from the master equation rewritten as a matrix problem of an associated 
quantum Hamiltonian H, the longest relaxation time was obtained from the real part of 
the lowest eigenvalue of H and steady-state observables were found by forming matrix 
elements with the right ground state of H. The results were compared with the ones 
of a Monte Carlo simulation. While the estimates of the location of the transition were 
in reasonable agreement between the two methods, it turned out that for the longest 
chains (L = 21 sites) accessible, finite-size estimates for the critical exponents could not 

|| He rather considered the process 2A — > 3^4, 3A — > A, which is believed to exhibit the same type of 
phase transition as the PCPD. 
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yet be reliably extrapolated towards the L — > oo limit [32]. As we shall see, even today 
the difficulties encountered in that first study are not yet completely overcome. 

At that time, the likely solution to these difficulties appeared to be the investigation 
of larger lattices by adapting White's density-matrix renormalization group (DMRG) 
[33,34] to non-equilibrium systems with non-hermitian matrices H (for a recent review, 
see [35]). Kaulke and Peschel [36] used the DMRG to study biased diffusion, which leads 
to a g-symmetric Hamiltonian H which is similar to a symmetric matrix. The feasibility 
of the Hamiltonian approach using the DMRG for truly non-equilibrium systems without 
detailed balance was demonstrated with E. Carlon and U. Schollwock, on a model in the 
DP class [37]. Equipped with this technique, we then returned to the investigation of a 
lattice model of the PCPD with an exclusion principle by density-matrix renormalization 
group methods [19]. This paper was followed by a long series of numerical and analytical 
studies [38-55] and released the still ongoing debate concerning the asymptotic critical 
behaviour of the PCPD at the transition. Until now a surprising variety of possible 
scenarios has been proposed, the main viewpoints being that the active-absorbing phase 
transition in the PCPD 

- should represent a novel universality class with a unique set of critical exponents [38, 
44,48,51], 

- may represent two different universality classes depending on the diffusion rate [39, 
52] and/or the number of space dimensions [53], 

- can be interpreted as a cyclically coupled DP and annihilation process [40] , 

- could be regarded as a marginally perturbed DP process with continuously varying 
critical exponents [47], 

- may have exponents depending continuously on the diffusion constant D [49,56], 

- may cross over to DP after very long time [50,54], and 

- is perhaps related to the problem of non-equilibrium wetting in 1+1 dimensions [55]. 

Given these widely different and partially contradicting conclusions, it is clear that 
considerably more work will be needed to understand the behaviour of this so simple- 
looking model. In this review, we shall give an introduction to the PCPD, review the 
research done on it and discuss these suggested scenarios in detail. We thereby hope to 
stimulate further research in order to finally understand the properties of this system. 

2. Phenomenological scaling properties 

It is well-known that systems in the directed percolation universality class are usually 
characterized by the so-called DP-conjecture [24, 30] which specifies certain conditions 
for DP critical behaviour. Similarly, experience suggests that the critical phenomenon 
observed in the PCPD is not restricted to a particular model, rather it is expected 
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to appear in a large variety of models which are thought to be characterized by the 
following phenomenological features: 

(i) Particles of a single species diffuse in a d- dimensional space. 

(ii) The phase transition is driven by two competing binary reactions for particle 
creation and removal, i.e., two particles have to come into contact in order to 
be eliminated or to produce offspring. 

(iii) There is a finite number of absorbing states of which at least one is reachable. 

(iv) There are no unusual features such as quenched randomness or unconventional 
symmetries. 

This means that the characteristic critical behaviour of the PCPD is also expected to 
exist in many other models, including e.g. generic reaction-diffusion processes of the 
form 

(a) 2A-+3A, 2A^$ 

(b) 2A^3A, 2A^A 

(c) 2A^AA, 2A^A 

(d) 2A -> 3A, 3A -> 

All these models differ significantly from both DP and PC processes in two respects: 
Firstly, processes in the DP class are unary processes, i.e. isolated particles can decay 
or trigger branching. The paradigmatic example is the contact process [4, 6] with the 
reactions A — > 2A, A — > 0, exclusion on the lattice, and possibly diffusion of single 
particles. In contrast, in the PCPD single diffusing particles cannot react unless they 
collide in pairs. Secondly, in a typical system in the PC class such as 2A — > 0, A — > 3A, 
the parity of the total number of particles is conserved and there is at least one unary 
reaction while this is not so in the PCPD. In particular, all PCPD-like models have an 
inactive (or absorbing) phase, where either the annihilation process 2A — > or else the 
coagulation process 2A — > A dominates. Therefore, in the inactive phase of the PCPD 
and related models the density p(t) of particles decays algebraically as [57] 

{t- d / 2 for d < 2 
t-Hnt forrf = 2 . (2) 
t- 1 for d > 2 

This behaviour is quite distinct from the time-dependence of the particle density for 
systems in the DP universality class, where the decay of p(t) is exponential in time. 

Studying the PCPD one has to distinguish between two fundamentally different 
types of models, namely, the unrestricted PCPD, where the number of particles per site 
is not constrained, and the restricted PCPD, where the particle number on any given site 
is effectively bounded.^ This restriction can be implemented either as a hard constraint 

% These two variants are often referred to as bosonic and fermionic models. This nomenclature, 
however, is misleading since e.g. models with a soft constraint are restricted but bosonic. 
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by an exclusion principle or as a soft constraint in form of a higher-order process for 
particle removal (such as 3 A — > 0), which prevents the particle density from diverging. 
The difference between the two types is most pronounced in the active phase, where in 
the unrestricted case the particle density diverges while it saturates at a finite value in 
restricted variants. Since phase transitions with a well-defined particle density in the 
active phase are perceived as more natural, most studies consider the restricted PCPD. 

In all recent studies of the restricted PCPD it is assumed that the critical behaviour 
at the transition can be described in terms of the same scaling laws as found in other 
absorbing phase transitions, although possibly with a different set of critical exponents 
and scaling functions. For example, in a Monte Carlo simulation one usually measures 
the temporal evolution of the particle density starting from a homogeneous initial state 
(usually a fully occupied lattice). In this case the particle density p is expected to obey 
the asymptotic scaling form 

p(t,L,e)=r s f^t^/L^t 1 ^) , (3) 

where L is the lateral system size, e = p — p c measures the distance from criticality, 
and f p is a scaling function with an appropriate asymptotic behaviour. Here and in 
what follows a scaling limit t — > oo, L — > oo and e — > is implied. The exponent 
8 — (3/v\\ describes the power-law decay of the density at criticality while z = v\\/v±_ is 
the dynamical critical exponent. Apart from metric factors the scaling function and the 
critical exponents are determined by the universality class of the phase transition and 
the boundary conditions. 

Alternatively one can start the process with an active seed, i.e. in the present case 
with a pair of particles at the origin. As usual in this type of simulation [23] , one studies 
the survival probability P s which scales as 

P s (t,L,e)=r s ' fp^LtV^et 1 ^) , (4) 

where 5' = 13'/ v\\. Similarly, the average number of particles averaged over all runs scales 

as 

N(t,L,e)=F faiLi 1 ';*; 1 '"*) , (5) 
and the mean-square spreading from the origin averaged over surviving runs is 

R 2 (t,L,e)=0 z ^(LfV^i/i) . (6) 

Here 5' = ft jv\\ is an independent exponent, although in the DP and the PC class the 
two exponents f3 and (3' are equal. + The so-called initial-slip exponent 77 is related to 
the other exponents by the generalized hyperscaling relation [3,58] 

77 = - - 8 - 8' . (7) 

z 

+ In DP the equality is caused by a time-reversal symmetry [23] while for the PC class it is so far only 
observed numerically 
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The definition of the survival probability in the PCPD is a little subtle. Clearly, the 
presence of a single particle as a criterion for survival is not sufficient since then P s (t) 
would tend to a constant. A natural definition, however, appears to be to consider a run 
as surviving as long the system has not yet reached one of the two absorbing states, i.e. 
there have to be at least two particles in the system. A stronger condition is to require 
the presence of at least one pair in the system as advocated in [47,49], see sections 4.3 
and 5.2. 

According to standard scaling theory the quartet of exponents (f3, (3', u±, v\\) 
characterizes the universality class of the transition. For DP and PC transitions, where 
(5 = f3', the remaining three independent exponents in one spatial dimension have the 
values 



The central open question about the PCPD is whether it belongs to one of the known 
classes or whether it represents a novel universality class with an independent set of 
critical exponents. 

3. Mean-field approaches 

Mean-field theories are a convenient tool in order to obtain a first orientation about 
the possible critical behaviour of a statistical system. To be specific, we shall formulate 
the PCPD as a lattice model, where each site may either be empty or be occupied by 
a single particle. The model evolves by random-sequential dynamics according to the 
rules [19] 



where < p < 1 and < D < 1 are the control parameters. We are therefore treating 
a restricted version of the PCPD and shall comment on the unrestricted version as 
discussed in [31] as necessary. 

All mean-field treatments must at some stage neglect some of the correlations 
present in a model. A systematic way of constructing a sequence of improving mean- 
field theories was devised by ben-Avraham and Kohler [60] and extending earlier work 
[61,62]. In their so-called (n, m)-approximation they work with clusters of n sites. The 
approximation is made by factorizing the probabilities for configurations which occupy 
more than n sites in terms of n-site probabilities but in a way that there is an overlap 
of m sites between adjacent clusters. More specifically, if 21, <B, <£. are such overlapping 
clusters, and if P(2l|23) = P(2l)/P(*8) is the conditional probability to find 21 for given 




{0.276486(8), 1.096854(4), 1.733847(6)} for DP [59] 
{0.922(5), 1.85(3), 3.22(6)} for PC [8] 



(8) 



AAQ, ®AA -> AAA 
AA — > 00 
A0 <-> 0A 



with rate (1 - p)(l - D)/2 
with rate p(l — D) 
with rate D , 



(9) 
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reaction 


AN, rate 


• • o — > • • • 
o • • — > • • • 

• • — > oo 


+1 \{l-p){l-D)p\l-p) 
+1 \(l-p){l-D)p\l-p) 
2 p(l - D)p 2 



Table 1. Rates in the site approximation of the PCPD. 



03, the approximation P(2UB(£) m P(2l<B)P(*B£|<B) is used. In this way, correlations up 
to n sites are taken into account. It is assumed throughout that the cluster probabilities 
P(2l) are spatially translation- invariant. 

For example, consider a five-site cluster in the state ABCDE. In the (3, 1)- 
and (3, 2)-approximations, respectively, the probability P(ABCDE) of this cluster is 
expressed as 

P^\ABCDE) » P{ABC)P{CDE\C) = P ( ABC ) P ( CDE ) 

P{C) 

P^(ABCDE) « P(ABC)P(BCDE\BC) » P(ABC) P ( BCD ) P ^DE\CD) ^ 

P(BC) 

P(ABC) P(BCD) P(CDE) 
~ P(BC)P(CD) ' 

Therefore, in the context of the (3, m)-approximation, the basic variables have to be 
constructed from the independent three-sites probabilities. 

It has been argued that the (n,n — ^-approximations are qualitatively the most 
reliable ones [60] and this conclusion has been generally accepted. In particular, the 
simplest of these schemes is the (1, 0)- approximation, which is also referred to as the 
site approximation (N = 1) or as simple mean-field. Similarly, the (2, l)-approximation 
is called the pair approximation (N = 2) and the (3, 2)- and (4, 3)-approximations are 
referred to as triplet (N = 3) and quartet approximations (N — 4), respectively. 

In what follows we shall denote an occupied site by • and an empty site by o. 
A three-site cluster AAfy will then be denoted as • • o and has the probability 
P(AA$) = P„ - For notational simplicity, we shall only consider the one-dimensional 
case in this section. 



3.1. Simple mean- field approximation for the particle density 

The simplest mean-field approach considers only the single-site probabilities. If we let 
pit) = P,(t) and also assume that it is space-independent, the rate equation for p(t) is 
easily derived. In table 1 we collect those reactions from (9) which change the mean 
particle-density, together with the change AiV. in the number of occupied sites N, by 
this reaction. 
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Adding the contributions to p and rescaling time t — > t(l — D) 1 one easily finds 

^ = {i-p) P (t?(i-p(t))-2 PP {ty (ii) 

which implies the following large-time behaviour of the mean density [19] 



^ + aexp(-t/r) 

p(t)~{ ^m-t-^ 

(3p- l)- 1 -r 1 



if p < 1/3 

if p = 1/3 (12) 
if p > 1/3 



where r = (1 — p)/(l — Sp) 2 and a is a constant which depends on the initial conditions. 
We see that the site approximation reproduces the intuitive expectation of a phase 
transition between an active phase, where p(t) p^ > for large times, and an 
absorbing phase, where p(t) vanishes in the t — > oo limit. We also see that throughout 
the absorbing phase, the approach towards the steady-state of vanishing particle density 
is algebraic, although mean-field theory is not capable of reproducing the correct decay 
p(t) ~ t^ 1 / 2 which one would obtain in one dimension if the fluctuation effects neglected 
here were taken into account. 

From the point of view of a continuum theory it is natural to introduce a coarse- 
grained particle density p(x,t) and to postulate a suitable rate equation. A reasonable 
rate equation should be given by* 

^p(x, t) = bp 2 (x, t) - cp 3 (x, t) + DV 2 p(x, t) , (13) 

and we expect for a restricted model with the rates (9) that b = 1 — 3p and c = 1 — p, 
see eq. (11). On the other hand, for unrestricted models without any constraint on the 
occupation of a space point x, one would have c = [31] so that the mean particle 
density p(x,t) diverges in the active phase b > 0. Hence the cubic term with c > can 
be interpreted to mean that a newly created particle requires some empty space to be 
put into. In later sections, we shall discuss variants of the PCPD with a soft constraint, 
generating effectively such a cubic term with c > 0. As we shall see this will lead to 
important consequences. 

At first sight, one might believe that eq. (11) contains with respect to (13) the 
further approximation that the space-dependence of p(x, t) \— > p(t) is neglected. Still, 
both eqs. (11) and (13) lead to the same long-time behaviour of the density. To see this, 
choose some spatial domain Q and consider the spatially averaged density 

m = ^J^ P (x,t), (i4) 

* We warn the reader that simple rate equations may not at all be adequate for the case of an 
unrestricted density, see [63] and section 6. 
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AN.. 


rate 






• • 


oo - 


-> • o •• 





- 1 


Duvwj p{\ — p) 




x2 


o • 


o» - 


-> o o •• 





+1 


Dv 3 /p(l - p) 




x2 


• • 


• • - 


-> • o o« 


- 2 


- 3 


p(l - D)u 3 /p 2 






• • 


• o - 


-> • o oo 


- 2 


- 2 


p(l - D)u 2 v/p 2 




x2 


o • 


• o - 


-toooo 


- 2 


- 1 


p(l - D)uv 2 /p 2 






• • 


o» - 




+1 


+2 


\{l-p){l-D)uv 2 /p{l- 


-P) 


x2 


• • 


oo - 


-> • • »0 


+1 


+1 


\{l - p){\ - D)uvw / p{\ 


-P) 


x2 



Table 2. Rates in the pair approximation of the ID PCPD. In the last column, a 
symmetry factor coming from parity-symmetry is included. 



where \ fl\ denotes the volume of fl and p(x, t) solves (13). Now, if c > and p(0) > and 
with the boundary condition Vp\on = 0, it can be shown rigorously [64] that there exists 
a constant b' so that b'/b > 0, and a positive constant c', which are both independent 
of \Q\ and of p(0) such that the inequalities 

bp(t) 2 - c'p(t) 3 < ^< b'p(t) 2 - cp{tf ; if 6 > 

(XL 

-\b'\p(t) 2 - c'p(t) 3 <^T< -NpW 2 - c P(tf ; if b < (15) 

hold true. In addition, one has in the quadratic mean the convergence p(x,t) — > p(t) 
for large times [64]. Therefore, the time-dependence of the averaged density p(t) can 
be bounded both from above and from below by a solution of the form (12). Hence 
we can conclude that the space-dependence of p(x,t) is not essential for the long-time 
behaviour of the mean density, rather the crucial approximation in mean-field theories 
is made in neglecting correlations between the states of different sites. 

Finally, a standard dimensional analysis of eq. (13), see [2], yields the mean-field 
critical exponents 

(3 MF = 1, v ± MF = l, v ll MF = 2. (16) 

As will be discussed in section 6, these exponents are expected to be valid in d > 2 
dimensions. In fact, high-precision simulations in 2D performed by Odor, Marques and 
Santos [45] confirm this prediction within numerical errors and up to logarithmic factors 
for various values of the diffusion rate. In one spatial dimension, however, fluctuation 
effects are expected to be relevant, leading to different values. A better description 
of these fluctuation effects requires to consider larger clusters, going beyond the site 
approximation. These approaches will be discussed in the following subsections. 



3.2. The pair approximation: Two transitions ? 

Treating the ID PCPD in the pair approximation, we assume again that the probabilities 
P(AB) are translation-independent and furthermore that the system is left/right- 
invariant, viz. P. (t) = P .(t). Then, because of P. = P. + P.. and P.. + 2P. + P 00 = 1, 
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there are two independent variables which may be chosen as the particle density 
p(t) = P.(t) and the pair density u(t) = P„(t). 

We now illustrate the standard methods [1,60] how to find the equations of motions 
for p{t) and u(t) in the pair approximation. Introduce the shorthand v = P, Q = p — u 
and w = Pqo = 1 — 2p + u. The reactions changing the number N, of occupied sites of 
the number N„ of occupied pairs are listed with their rates in table 2, grouped into the 
separate contributions of the three reactions. In the last group, we only need to take 
into account those reactions 2A — > 3A which modify the particle configuaration on a 
given site. Adding the contributions to p and ii, their equations of motion are easily 
derived [19]: 

p(t) = - 2(1 -D)pu(t) + (I -D)(l-p)(p(t)-u(t))^ (17) 

u( t) = - (i - DH^wyw - u>w - - ^t) 2 ) 



+ (l-D)(l-p) 



P(t) p(t)(l-p(t)) 
(p(t) - u(t))(l - u(t))u(t) 



p(t)(l-p(t)) v ' 

The pair approximation is the lowest-order cluster approximation which allows the 
effects of diffusion to be treated explicitly. In particular, in the D — > 1 limit, the 
site approximation eq. (11) is recovered for p(t), while the pair density u(t) = p(t) 2 . 

The critical line is again at p c (D) — | for D > y, but if < D < y, one has 
p c {D) = |(1 + 3D) /(I — D). In addition, the scaling of both p(t) and of u(t) - and 
similarly of their steady-state solutions p^ and - does depend on whether D > | 
where u ~ p 2 or whether D < y where u ~ p. The steady-state particle density 
Poo ~ Pc(D) —p for all D, and if D ^ j one has u(t) ~ t^ 1 along the critical line. These 
results might suggest the existence of two distinct universality classes along the critical 
line, but since the pair approximation is the simplest cluster mean-field theory which 
permits to study the effect of particle diffusion at all, this conjecture certainly needs 
further verification. 

In figure 2 the phase diagram of the ID PCPD is shown. The pair approximation 
(N = 2) certainly is an improvement as compared to simple mean- field (N = 1) but is 
still quite far from the exact location of p c {D) as obtained from DMRG calculations or 
Monte-Carlo simulations which will be described in section 4. 



3.3. Higher cluster approximations 

In the pair approximation of the PCPD, the particle and pair densities in the steady- 
state are related by Uoo = poo(l — 3p)/(l — p). This means that the particle and pair 
densities have to vanish simultaneously which is certainly not true for the simple PCP 
without diffusion (D = 0), see section 1. This failure of the pair approximation can 
be overcome within the triplet approximation (N = 3). Consider first the PCP where 




Figure 2. Phase diagram of the one-dimensional restricted PCPD with the rates as 
defined in eq. (9), after [52]. The figure shows the boundary between the active and 
the absorbing phases according to cluster mean- field theory, ranging from site (N = i) 
to quartet approximation (N — 4). The filled dots give values of p c (D) as found from 
the DMRG [19] and the open triangles those found in Monte Carlo simulations [65]. 
The grey diamond gives the critical point p c (0) of the PCP [21]. 



D = [22,66]. Again assuming left/right symmetry, the independent variables may 
be chosen as P..., P,, , P» 00 and P. .. Now the critical point occurs at p c (0) — 0.128, 
where the steady-state pair density Uoo ~ p c (0) — p vanishes linearly, but the critical 
particle density ~ 0.23 remains finite. These results have been extended up to the 
sextet approximation (N = 6) [67] and further to N = 12 [68] and remain qualitatively 
valid. 

It is then natural to extend this study to the whole phase diagram. This was done 
in [52], including the quartet approximation (N = 4). As can be seen in figure 2 the kink 
in p c {D) observed in the pair approximation is absent for larger clusters. In addition, 
along the critical line, the steady-state particle density ~ p c {D) — p for all values 
of D > 0, whereas ~ (p c (D) — p) 2 and this holds for both A = 3 and N = 4. 
This result strongly suggests that there is only one universality class along the critical 
line [52]. This finding is in agreement with simulational results, see section 4. A similar 
result had already been obtained before in the 2D process 2A — > 0, 2A — > AA with 
single-particle diffusion [45]. 

To summarize, the apparent presence of two universality classes for D > should 
be an artifact of the pair approximation - at least in one dimension - and it has already 
been observed earlier that in some cases, cluster mean-field theories may indeed yield 
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qualitatively incorrect phase diagrams [60]. Generically, in the PCPD some reactions 
involve three neighbouring sites, and furthermore the explicit consideration of at least 
three neighbouring sites is required in order to become sensitive to the large number of 
absorbing states. Therefore it appears natural that cluster mean-field approximations 
which already neglect some of the correlations present in the reaction terms themselves 
are unlikely to be reliable (in fact, the site approximation may be too simplistic to yield 
anything but a single transition). Rather, the pair approximation is good enough to 
reveal the presence of a PCP transition for D = which is distinct from the one which 
occurs for D > 0, but since three-sites correlations are not treated correctly, the PCP 
transition might have been stabilized artificially even for some values D > 0. Similar 
considerations should apply to any interacting particle system where the interactions 
occur through the contact of at least two particles (see also section 6). 

4. Numerical studies of the critical behaviour in one dimension 

4-1. Density-matrix renormalization group study 

Since the first non-perturbative study of the PCPD [19] was made using the 
reformulation of the problem through the quantum hamiltonian formalism, we shall 
begin with a brief description of that method, see [2, 69, 70] for recent reviews. The 
starting point is provided by the master equation, cast into the form 

d t \P(t)) = -H\P(t)) , (19) 

where \P(t)) = J2{a} P({ a }'^)\{°'}) * s the state vector of the system at time t, \{cr}) 
is the state vector of a particle configuration {a} on the lattice, and P({a};t) is the 
time-dependent probability of that configuration. For a chain of L sites with a maximal 
occupancy of one particle per site the quantum Hamiltonian H is a 2 L x 2 L matrix with 
elements 

({a}\H\{r}) = -w(t - a) + S a , T £ w(r - v) , (20) 

M 

where w(t — > a) is the transition rates from the configuration {r} to the configuration 
{a}. Since the elements of the columns of H add up to zero, H describes indeed a 
stochastic process. From the condition of probability conservation P({a}; t) — 1 it 
follows that (s\H = 0, where (s\ = X^o-} i s a eigenvector of H with eigenvalue 0. 
Since the particle reactions are generically irreversible, H is in general non-hermitian and 
therefore it has distinct right and left eigenvectors, which we denote by |0 r ), |l r ), . . . , \n r ) 
and (0^1, (1/|, . . . , (ni\, respectively with the corresponding 'energy levels' E Q , Ei, . . . , E n 
ordered according to = Eq < Re E\ < • • - Re E n . Clearly, (n/ \m r ) = if E n ^ E m 
and we normalize states such that (ni\n r ) = 1. We also set (s\ = (0i\. 
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Given the rates (9) of the PCPD, the associated quantum Hamiltonian H on a 
chain with L sites is easily constructed. The stationary states of the model correspond 
to right eigenvectors of H with vanishing eigenvalue. For the PCPD with D > these 
are, for both periodic and free boundary conditions [19], the two absorbing states 

|0 r ) := |o o o • • • o) 

|l r ) : = 1 (|« o o • • -o) + |o • o • • -o) H h|ooo--.«)) . (21) 

Consequently, the first two eigenvalues E and E 1 vanish so that the inverse relaxation 
time towards the absorbing state is given by the smallest 'energy gap' T = T(p, D; L) : = 
Re E 2 . Finally, formally solving the master equation (19), the average value (£>}(t) of 
an observable O is given by 

(£>)(*) = (*\D\P(t)) =J2( s \°\ n r)( n i \m)e~ Ent , (22) 

n 

where we tacitly assumed the completeness (and biorthogonality) of the right and 
left eigenvector systems. Although the quantum hamiltonian formalism was originally 
introduced with the aim of relating non-equilibrium systems on infinitely long chains to 
integrable quantum systems, we are interested here in the numerical analysis of systems 
on finite chains with L sites. Then the long-time behaviour is given by the smallest 
inverse relaxation time r = £y = with the following finite-size behaviour 



exp(-L/£j 

L- 2 



if p < p c (D) 

iip = p c (D) (23) 
if p > p c {D) 



where z = ^||/V_i_ and £j_ ~ \p — p c (D)\~ u± is the spatial correlation length. For large 
lattices and p > p c (D), one expects that the lowest energies should satisfy a massless 
dispersion relation of the form E(k) ~ k z , and since the lowest momenta should scale 
as /c m in ~ L -1 , the phenomenological scaling (23) is recovered. 

Having found T(p,D;L) on a sequence of finite lattices, one proceeds as usual, 
see [71] for a collection of reviews. The critical parameters (here p c (D) and z) are 
extracted by forming the logarithmic derivative 

Yl ^ V >- ln[(L + l)/(L-l)] • [M) 

However, since the entire absorbing phase is critical and therefore lim^oo F(p, D; L) =0 
the habitual method of looking for the intersection of two curves Y L (p, D) and Y L i(jp, D) 
cannot be applied [19]. Rather, the critical point p c {D) must be found by fixing d and 
L and finding the value of p = Pl{D) which maximizes Yl(j>,D) (see Appendix B for 
further details). Then the sequence of estimates has to be extrapolated to p c (D) = 
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lim 



l^ooPl(D). Finally, for L large, one expects 



if p < p c (D) 
ifp = p c {D) 
if p > p c (D) 



Y L (p,d) 



z 



(25) 



2 



and one finally obtains z = — Hindoo Y L (p c (D), D) [19]. In addition an estimate of the 
order parameter exponent (3/v± can be obtained from the steady-state density profile 
at p = p c {D) 



where n(£) is a particle number operator at site i and / is a scaling function. The 
double degeneracy of the ground state was lifted by adding a particle creation process 
— > A with a small rate p'. For estimates of the bulk exponent, one may set I = L/2 
and p L (L/2) ~ L"^ as expected. On the other hand, if one sets £ = 1, in the limit 
of p' — > one has the scaling form Pl(1) ~ LiPi-P)/^ from which the surface critical 
exponent /3i/v± can be determined [37] .(j 

The lowest eigenvalue T and the corresponding eigenvector are found from the non- 
symmetric DMRG algorithm [19,37] on an open chain. In comparison with Monte-Carlo 
simulations, the DMRG method as described here directly yields the steady-state with 
a high numerical precision, does not suffer from a critical slowing-down in the vicinity of 
the critical point and for non-disordered systems does not require the use of a random- 
number generator. On the other hand, the number of lattice sites which can be treated 
is relatively restricted. For example, in [19] data for Tl with L < 30 sites and for 
with L < 48 sites could be obtained. Larger lattices with up to L = 60 sites have been 
obtained more recently [44,54]. In the absorbing phase, however, the finite-size sequences 
converge very well and the accuracy of the estimated limits can be considerably enhanced 
through the bst extrapolation algorithm [74,75]. The expected dynamical and order- 
parameter exponents z = 2 and (3/v± = 1 are recovered to high accuracy. On the 
other hand, along the critical line, the required finite-size extrapolations are affected by 
some ill-understood correction terms which preclude the use of sequence extrapolation 
algorithms. The final results for p c {D) and the exponents z and /3/v± are collected in 
table 3. Since results for several values of D were obtained, a discussion of systematic 
effects on the estimates becomes possible. 

While the values of /3/u± are almost independent of D, the results for z show 
a considerable variation which was argued in [19] to come from a subtle finite-size 
correction. Taking mean values and comparing with the values of the other universality 
classes available at the time, namely directed percolation (DP) and the parity-conserved 

jj For models with a non-degenerate ground state the use of a surface rate may be avoided by estimating 
the profile from N(£) = (l;|n(£)|l r ). At criticality, the finite-size scaling forms N(L/2) <~ L~ fi l Vl - and 
N(l) - L-/V^ are expected [37,72,73]. 




(26) 
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(PC) universality classes, the values of z and f3/v± of [19] are far from those of DP but 
intriguingly close to those of PC, see table 3. Because of this surprising coincidence 
the possibility that the active-inactive transition of the PCPD might belong to the PC 
universality class emerged, in spite of the fact that there is no obvious counterpart for 
the parity conservation law in the PCPD, i.e. there is no symmetry separating the 
dynamics into two sectors [38]. In fact, in both cases the inactive phase is characterized 
by an algebraic decay of the density as t^ 1 / 2 . Moreover, both models have two absorbing 
states. Still, the unambiguous identification of a steady-state universality class requires 
the determination of four independent exponents, see [1-3], rather than the two obtained 
in [19]. 

Indeed, small-scale Monte-Carlo simulations [38,39] seemed to confirm the DMGR 
results for p c (D) and, for D not too large, also the exponents z and (3/v±. However, 
either the exponent bound f3 < 0.67 [38] or the estimate (3 ~ 0.58 [39] was significantly 
smaller than the expected PC value 0.92. In addition, the effective exponents 5' and rj in 
seed simulations were shown to differ significantly from the expected PC values. These 
results rather suggested that the PCPD might be a novel nonequilibrium universality 
class, distinct from that of PC as well as DP processes. 

4-2. High-precision Monte- Carlo simulations 

The first simulational study of the PCPD with sequential dynamics (the rates (9) were 
chosen for maximal efficiency) yielded the critical-point estimates [77] 

p c (0.1) = 0.1105(5) , p c (0.5) = 0.153(1) (27) 

which are in good agreement with later results shown in figure 2 and table 3. We point 
out that several inequivalent definitions of the critical threshold p c (D) are in use and 
only the definitions of the articles [19,38,47,54] and of figure 5 - which gives the results 
of [65] - are consistent with each other and eq. (27). A systematic high-precision Monte- 
Carlo simulation was performed by Odor [39]. In order to implement the model on a 
parallel computer, he investigated a variant of the PCPD with synchronous updates on 
three sublattices. Depending on the value of the diffusion rate he found continuously 
varying exponents 5, (3 and rj (see table 3), concluding that the process does not belong 
to the PC universality class. Instead the value of (3 seemed to jump from (3 0.58 
for 0.05 < D < 0.2 to smaller values at D > 0.5. This observation together with 
the aforementioned pair mean-field approximation led Odor to the conjecture that the 
PCPD might belong to two different universality classes for low and high diffusion 
rate. Although it seems to be surprising how the choice of a rate can switch between 
two universality classes we note that a similar scenario has been observed in a certain 
two-species particle model [78], although the specific mechanism in these models is not 
obvious in the PCPD. We shall return to this problem in section 6. 
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Reference 


D 


Pc 


5 


^11 




0/"± 


5' 


V 


Carlon et al. [19] 


0.1 


0.111(2) 


— 


1.87(3) 


— 


0.50(3) 


— 


— 




0.15 


0.116(2) 




1.84(3) 




0.49(3) 








0.2 


0.121(3) 


— 


1.83(3) 


— 


0.49(3) 


— 


— 




0.35 


0.138(1) 




1.72(3) 


— 


0.47(3) 


— 


— 




0.5 


0.154(1) 




1.70(3) 


— 


0.48(3) 


— 


— 




0.8 


0.205(1) 




1.60(5) 




0.51(3) 






Hinrichsen [38] 


0.1 


0.1112(1) 


0.25(2) 


1.83(5) 


< 0.67 


0.50(3) 


0.13(2) 


0.13(2) 


Odor [39] 


0.05 


0.25078 


0.273(2) 


— 


0.57(2) 


— 


0.004(6) 


0.10(2) 




0.1 


0.24889 


0.275(4) 


— 


0.58(1) 




— 


— 




0.2 


0.24802 


0.268(2) 




0.58(1) 




0.004(6) 


0.14(1) 




0.5 


0.27955 


0.21(1) 




0.40(2) 




0.008(9) 


0.23(2) 




0.9 


0.4324 


0.20(1) 




0.39(2) 




0.01(1) 


0.48(1) 


Odor [41] 


0.2 


0.17975(8) 


0.263(9) 




0.57(1) 










0.4 


0.2647(1) 


0.268(8) 




0.58(1) 










0.7 


0.3528(2) 


0.275(8) 




0.57(1) 








Park et al. [42] 


* 


0.0895(2) 


0.236(10) 


1.80(2) 


0.50(5) 




w 0.1 


w 0.2 


Noh & Park [47] 


0.1 


0.1112(1) 


0.27(4) 


1.8(2) 


0.65(12) 


0.50(5) 


0.09(2) 


0.18(5) 


Park & Kim [48] A 




0.03081(4) 


0.241(5) 


1.80(10) 


0.519(24) 


— 


0.11(3) 


0.15(3) 


Park & Kim [48] B 




0.28735(5) 


0.242(5) 


1.78(5) 


0.496(22) 




0.13(3) 


0.15(3) 


Dickman 


0.1 


0.10648(3) 


0.249(5) 


2.04(4) 


0.546(6) 


0.503(6) 






OUT r a r\ 1 

& Menezes [49] 


0.5 


0.12045(3) 


0.236(3) 


1.86(2) 


0.468(2) 


0.430(2) 








0.85 


0.13003(1) 


0.234(5) 


1.77(2) 


0.454(2) 


0.412(2) 






Hinrichsen [50] 


1 


0.29998(4) 


0.22(1) 


1.78(5) 






0.13(1) 


0.22(3) 


Kockclkoren 


* 




0.200(5) 


1.70(5) 


0.37(2) 


— 


— 


— 


& Chate [51] 


















Odor [52,56] 


0.05 


0.10439(1) 


0.216(9) 


2.0(2) 


0.411(10) 


0.53(7) 


— 


— 




0.1 


0.10688(1) 


0.206(7) 


1.95(1) 


0.407(7) 


0.49(2) 


— 


— 




0.2 


0.11218(1) 


0.217(8) 


1.95(1) 


0.402(8) 


0.46(3) 








0.5 


0.13353(1) 


0.206(7) 


1.84(1) 


0.414(16) 


0.41(2) 








0.7 


0.15745(1) 


0.214(5) 


1.75(1) 


0.39(1) 


0.38(2) 






Barkema 


0.1 


0.11105(5) 


0.17 












& Carlon [54] 


0.2 




0.17 


1.70(3) 




0.28(4) 








0.5 


0.15245(5) 


0.17(1) 






0.27(4) 








0.9 


0.2335(5) 


0.17 


1.61(3) 










Related models: 


Hinrichsen [40] 


1 


0.6929(1) 


0.21(1) 


1.75(10) 


0.38(6) 




0.15(1) 


0.21(2) 


Odor [46] 


* 


0.3253(1) 


0.19(1) 


1.81(2) 


0.37(2) 








Hinrichsen [76] 


* 




0.21(2) 


1.75(2) 










Known universality classes: 


DP [59] 






0.1595 


1.5807 


0.2765 


0.2521 


0.1595 


0.3137 


PC [8] 






0.286(2) 


1.74(1) 


0.922(5) 


0.497(5) 





0.286(2) 



Table 3. Summary of estimates for the critical exponents of the II? PCPD. A star 
indicates that the control parameter is effectively fixed or defined in a different way. 
Several variants of models with differently defined critical thresholds p c {D) were used. 
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As pointed out in [38] the unusual critical behaviour of the PCPD should be seen in 
a large variety of binary reaction-diffusion processes with absorbing states, in particular 
in the coagulation-fission process 2A — > 3A, 2A — > A. Simulating the coagulation-fission 
process with branching sidewards 

AA®, ®AA -> AAA with rate (1 - p) (1 - D)/2 

AA -> 0A/A0 with rate p(l - D) /2 (28) 

A0 <-> 0A with rate D 

as well as a symmetric variant with offspring production in the middle of two particles 

A$A -> AAA with rate A(l - p){l - D) (29) 

Odor [41] could confirm this expectation, provided that the diffusion rate is high enough. 
The observed exponents seemed to be in agreement with those of the PCPD for small 
diffusion rates, leading Odor to the conjecture that the coagulation-fission process may 
belong to a single universality class. 

The special case, where the rates for coagulation and diffusion are equal, i.e. 
D = p(l — D)/2, can be solved exactly for any value of A [43], using the well-known 
method of interparticle distribution functions [79]. As pointed out in [43], in this case 
the dynamics is special in so far as the fission process can only fill voids of an invariant 
skeleton of coagulating random-walks. Unfortunately, in this case the active phase is 
inaccessible so that the exact solution does not provide any information about the critical 
exponents of a possible universality class of PCPD transitions, rather the system always 
crosses over to a decay of the form pit) ~ t -1 / 2 . Nevertheless the solution reveals that 
the finite-size scaling of the relaxation times does not depend on the fission rate and 
thereby allows to test a recent extension to non-equilibrium phase transitions [44] of the 
Privman-Fisher finite-size scaling forms [80]. 

Concerning the possibility of two universality classes or even continuously varying 
exponents, Park and Kim [48] showed by numerical simulations that in models, where 
the rates of diffusion and annihilation are tuned in such a way that the process without 
particle production would become exactly solvable, one obtains a well-defined set of 
critical exponents [48] (see table 3). In this sense, they argue, the transition in the 
PCPD may be considered as universal. Contrarily Dickman and Menezes [49] claim 
that the transition in the PCPD is non- universal. Their argument will be reviewed in 
the next subsection. 

As a third possibility one of us (H.H.) suggested an extremely slow crossover to 
DP [50] . Simulating a cellular automation of the PCPD which was originally introduced 
by Grassberger [65], it was argued that the scaling regime is not reached, even after six 
decades in time. Instead the effective exponents display a slow drift as time proceeds. 
Remarkably, all exponents seem to drift in the direction of DP values so that one cannot 
rule out that the critical behaviour of the 1+1-dimensional PCPD may eventually tend 
to DP. According to these arguments the DP process manifests itself in the dynamics of 
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Figure 3. Left panel: Decay of the particle density in the inactive phase for different 
values of p in comparision with the critical PCPD (bold line) and pure coagulation 
(dotted line). Right panel: Sketch of the sudden decay of the particle density in the 
inactive phase (see text). After [50]. 



pairs of particles, while the diffusing background of solitary particles becomes irrelevant 
as t — > oo (obviously, this mechanism should work exclusively in one spatial dimension 
while in higher dimensions the diffusive background may be increasingly relevant leading 
to an effective mean- field behaviour). Recently Barkema and Carlon [54] presented 
extrapolation results, arriving at a similar conclusion (see section 4.4). Contrarily, a 
recent simulation of the PCPD with non-restricted occupancy per site performed by 
Kockelkoren and Chate (see section 4.5) does not show such a drift, supporting the 
hypothesis of a single universality class distinct from previously known classes. 

The DP hypothesis poses a conceptual problem, namely, how does the algebraic 
decay pit) ~ t~ 1 / 2 in the inactive phase comply with DP, where the decay is exponential? 
In this context it is interesting to monitor the crossover of the particle density in a slightly 
subcritical PCPD. As shown in figure 3, there are three different temporal regimes. 
In the first regime (region I) the particle density decays very slowly as in a critical 
PCPD. After a characteristic time, however, the density decays rapidly (region II) until 
it crosses over to an algebraic decay t -1 / 2 (region III). Remarkably, the amplitude of this 
asymptotic power law does not depend on p — p c . Therefore, the breakdown in region II 
becomes more and more pronounced as the critical threshold is approached and might 
tend to a quasi-exponential decay in the limit p — > p c . 

4-3. Universal moments and amplitudes 

Since the numerical values of the available exponents seem to vary between different 
universality classes, alternative diagnostic tools might be of value. Following [20], let us 
consider the order-parameter moments 



m n := (p n ) 



(30) 
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d 


model 


7714/7712 




to 3 /(toito 2 ) 


m 2 /m\ 


Ki/Kl 


Rcf. 


1 


DP 


1.554(2) 


1.526(3) 


1.301(3) 


1.1736(2) 


-0.505(3) 


[20] 




PCP 


1.558(2) 


1.529(3) 


1.303(3) 


1.1738(2) 


-0.493(3) 


[20] 




PC 








1.3340(4) 




[49] 


2 


DP 


2.093(8) 


2.080(1) 


1.569(1) 


1.3257(5) 


-0.088(4) 


[20] 




PCP 


2.07(1) 


2.067(9) 


1.56(1) 


1.323(3) 




[21] 



Table 4. Some universal moment and cumulant ratios at criticality for the directed 
percolation universality class (DP), the pair-contact process (PCP) and the parity- 
conserving class (PC) for an infinite lattice in d = 1,2 dimensions. 

From the scaling arguments reviewed in section 2 it is known that ratios such as m 4 /m 2 , 
m 3 /ml, m3/(m 1 m 2 ) or m 2 jvn\ should be universal. In fact, moment ratios are just 
special values of the scaling functions and hence their numerical values may be used, 
like the critical exponents, to identify the universality class of a given model. Further 
one may consider the cumulants [20] 

K 2 := m 2 — m 1 , K 4 := m 4 — Am^mi — 3m 2 + 12m 2 m 1 — Qm 1 (31) 

whose ratios K±jK\ and K 2 /ml = m 2 /m 2 — 1 should be universal as well. 

For future reference, we reproduce in table 4 several universal moment ratios 
as obtained for the simple contact process (DP), for the pair-contact process (PCP) 
without diffusion, and for a branching-annihilating random walk belonging to the parity- 
conserving class (PC). The numerical estimates were obtained at the critical point of 
finite-size systems by averaging over surviving quasi-stationary runs and extrapolating 
L — > 00. Note that for the contact process, the order parameter is the particle density, 
whereas the pair-contact process with D — 0, the order parameter is the pair density. 
While it had been known before that the static critical exponents of the DP and the 
PCP agree [17] the agreement of the moment ratios between the two models provides 
additional evidence in favour of the conjecture that the steady-state transitions of the 
DP and the PCP are in the same universality class. On the other hand, the result for 
the branching-annihilating random walk (BAW) [8] which is in the PC class is clearly 
different, as expected. We are not aware of analogous data for any other interacting 
particle model but the building-up of a collection of data of these amplitudes in several 
universality classes should be helpful. 

Having established that the method of universal moments is capable to distinguish 
between different universality classes, Dickman and Menezes [49] tried to apply the 
moment ratio method to the ID PCPD and considered the universal ratio m : = 
(p 2 )/(p) 2 = m 2 /m 2 . Although several critical exponents could be determined from 
finite-size scaling, surprisingly they found for the moment ratio a logarithmic increase 
with the system size asm~ (lnL)^ with ijj > 1, see figure 4. No such behaviour is seen 
for the BAW. If this phenomenon persisted in the limit L — > 00, it would collide with 
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Figure 4. Left panel: Probability distribution of the number of particles n at 
criticality, for D = 0.1 and +: L = 80; x: L = 160; □: L = 320. The inset 
shows m = m(lnL) of the PCPD and also for the BAW. Right panel: scaling of the 
reduced probability distribution in the reactive sector and with the same symbols as 
before. The inset shows m = m(lnL) in the reactive sector of the PCPD, with ■: 
D = 0.1; +: D = 0.5; x: D = 0.85. After [49]. 
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the scaling forms reviewed in section 2 from which one would deduce a scaling form 

P(p;L) = (p)V(p/(p}) , (32) 

where V is a normalized scaling function. Indeed, as can be seen in figure 4, the particle 
distribution function does not scale (a similar result also holds for the pair distribution 
function). Instead, the most probable particle number is always 2 and the distribution 
shows a tail which grows with L and generates the critical behaviour of the model. 

In order to account for this, Dickman and Menezes [49] argue that the PCPD 
should be characterized by the simultaneous dynamics of isolated particles and colliding 
pairs (as is indeed rendered plausible by illustrations of the dynamics as in figure 1) 
and therefore propose to separate the model's states into a 'purely diffusive' sector 
and a 'reactive sector' which must contain at least a pair of particles. The idea is to 
restrict the analysis to the reactive sector by excluding all states without pairs when 
taking averages, but without modifying the dynamical rules (this was also proposed 
independently in [47]). The result is shown in the right panel of figure 4, where a 
scaling behaviour in the reduced variables p* = pj (p) and P* = (p) P is obtained. Now, 
the form of P* becomes quite similar to the one of the non-diffusive PCP. In addition, 
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the ratio m also becomes almost constant and with a value close to the estimates for 
DP, see figure 4. Still, as already noted in [49], neither the scaling collapse nor the 
L-independence of m is completely perfect and the origin of these deviations remains to 
be understood.! f Although off-critical simulations lead to a clear scaling collapse and 
a DP-like exponent u± = 1.10(1) independently of D [49], other exponents maintain 
a weak D-dependence. When taking these results at face value, it implies a non- 
universality along the critical line and both the DP and the PC universality classes 
are excluded for the ID PCPD [49] . However, one might wonder whether the estimates 
of the exponents could still be affected by some non-resolved correction-to-scaling term 
which would make the dependence on D only apparent. In any case, this study points 
to the large corrections to scaling which are present in the ID PCPD and which render 
the extraction of universal parameters very difficult. 

While the above universal moment ratios referred to a spatially infinite system, we 
now consider the case of a finite lattice of linear size L, e.g. a hypercube. For equilibrium 
systems, it is well-known that the finite-size scaling of the correlation length leads to 
universal amplitudes, as reviewed in [81]. For non-equilibrium phase transitions, it has 
been argued that on a finite spatial lattice, the temporal correlation length ^ and the 
spatial correlation length should scale as [44] 

q^L-DoRi {C^p-pjL 1 ^) 

^ i = L- 1 S i {C 1 (p-p e )L 1 ^) (33) 

where (for given boundary conditions) Ri and Si are universal scaling functions, the 
index % refers to the observable Dj of which the correlators are studied and the entire 
non- universality can be absorbed into the non-universal metric factors C\ and D . 
In particular, the existence of the above scaling forms can be derived for systems in 
the directed percolation universality class [44]. For the convenience of the reader, 
the argument will be reproduced in appendix C. Consequently, the finite-size scaling 
amplitudes of the spatial correlation lengths ^J 1 are universal, as are ratios of temporal 
correlation lengths £|jjV£jj"j = ^i/Ej which in turn are easily calculated from the 
spectrum of the quantum hamiltonian H, as discussed above in section 4.1. Furthermore, 
from the scaling forms derived in appendix C it follows that the universality of the critical 
finite-size scaling amplitudes of £±j is equivalent to the universality [44] of the moments 

L L) d d r ± j™d n rlG(r^ n] L/U) 
/n(L) d r -LJo dr||G(r ± ,r||;L/^j.) 

of the pair connectivity G(r±,r\\) which is the probability that the sites (0,0) and 
(r±,r\\) are connected by a direct path. Here Q(L) is a d- dimensional hypercube of 
linear extent L and S is a universal function the argument of which does not contain 
any non-universal metric factor. These universality statements have been confirmed in 
several exactly solvable models [43,44]. In figure 5 we illustrate the universality of the 

ff Restricting to states with at least 2 pairs does not reduce further the apparent variation of m with L. 
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Figure 5. Universal ratio of finite-size scaling amplitudes of the relaxation times along 
the critical line in the ID PCPD, for free boundary conditions. The values of p c (D) 
are from [65]. After [44]. 

amplitude ratio R = E 3 /E 2 along the critical line in the ID PCPD, as obtained from 
the DMRG. As observed before for the exponents, there appear to be some finite-size 
corrections which are sufficiently irregular such that the powerful sequence-extrapolation 
methods [75, 82] cannot be used. Still, inspection shows that for the smaller values of 
D (i.e. up to D < 0.5) the data seem to converge towards a limit Roo xs 2 and its 
independence of D would confirm universality. On the other hand, for values of D > 0.8, 
there appears to arise a large transient region and only for values of L larger than ps 50 
sites a convergence might slowly set in. This is a cautionary example illustrating again 
that the truly asymptotic regime in the PCPD might set in very late and sometimes 
could be beyond the reach of currently available numerical methods. 

At present, no finite-size data for the amplitude ratios (33) or (34) seem to be 
available for free boundary conditions in any nonequilibrium universality class. In this 
context we recall the well-known observation from critical equilibrium systems that 
universal amplitude ratios tend to vary considerably more between different universality 
classes than usually do critical exponents, see [81]. Finding reference values for these 
amplitudes in distinct non-equilibrium systems should be helpful. 

4-4- Extrapolation techniques 

In numerical simulations of the PCPD there are two sources of uncertainty which cannot 
be separated easily. On the one hand, the critical PCPD shows unusually strong intrinsic 
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corrections to scaling, on the other hand a possible error in the estimation of the critical 
point leads to additional off- critical deviations with an unknown functional form. Most 
numerical approaches assume that the first source of uncertainty can be neglected after 
sufficiently long simulation time. The critical point is then estimated by postulating 
power-law behaviour over the last few decades in time, indicated by a straight line in 
a log-log plot or saturating effective exponents. In the PCPD, however, the main risk 
lies in the possibility that the intrinsic corrections to scaling at criticality may still be 
present even after 10 6 . . . 10 8 time steps. Thus, looking for a power-law over the last few 
decades of time, one may be tempted to 'compensate' the intrinsic scaling corrections 
by off-critical deviations in opposite direction, leading to a slightly biased estimate of 
the critical point and thereby to considerable systematic errors in the estimated critical 
exponents. 

In order to solve this problem the two sources of corrections - intrinsic scaling 
corrections and off-critical deviations - have to be separated. One may try to do this 
by assuming a certain functional form for one of them and performing appropriate fits. 
The first attempt in this direction was made by Odor [52], who postulates logarithmic 
corrections to scaling at criticality. For example, in the case of a temporally decaying 
particle density, where we expect an asymptotic power law p(t) ~ t~ s , he assumes the 
intrinsic scaling correction at criticality to be of the form 



where a and b are fit parameters. Using this technique Odor obtained 5 = 0.21(1) for 
various values of the diffusion rate between 0.05 and 0.7, favouring the scenario of a 
single independent universality class. Odor also estimated the critical exponents in a 
two-dimensional PCPD, confirming the mean-field prediction in d > 2 dimensions. 

Another extrapolation method was tried by Barkema and Carlon [54] in analysing 
high-quality data of a multispin Monte Carlo simulation. Unlike Odor, they postulate 
algebraic scaling corrections of the form 



where a and 7 have to be fitted appropriately. Although one would naively expect the 
exponent of the correction 7 to be equal to 1, Barkema and Carlon present numerical 
evidence suggesting that 7 = 5 = /3/v±. Based on this conjecture they propose to 
plot local slopes versus the density of particles. Using this representation the curve of 
local effective exponents should intersect with the vertical axis at p = as a straight 
line, as sketched in figure 6. Applying this extrapolation method Barkema and Carlon 
argue that the critical exponents tend to those of DP, in agreement with the conclusions 
offered in [50] (they only quote the value of 5 averaged over several values of D which 
is the one listed in table 3). 




(36) 
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Figure 6. Extrapolation technique suggested by Barkema and Carlon (after [54]). 
Assuming algebraic corrections of the form (36) with 7 = S the curve of the effective 
exponent 6 c s(t) should intersect with the vertical axis as a straight line (lower bold 
line). Ordinary extrapolations such as in Ref. [51], assuming that the scaling regime is 
actually reached within the simulation time, would correspond to the upper bold line, 
which becomes horizontal as p(t) — > 0. 



4-5. Restricted PCPD with a soft constraint 

All variants of the PCPD discussed so far are restricted by hard-core exclusion, i.e. each 
site can be either empty or occupied by a single particle. Although these models are 
easy to implement numerically, they involve three-site interactions with some ambiguity 
concerning the spatial arrangement. In order to circumvent this problem, Kockelkoren 
and Chate [51] proposed a surprisingly simple 'bosonic' model in which a soft constraint 
prevents the particle density in the active phase from diverging. The scaling arguments 
of the previous sections are expected to go through and models with a soft constraint 
are expected to exhibit the same type of critical behaviour as the ordinary PCPD with 
hard-core exclusion, see also section 6. 

Their model is defined as follows. Particles of a single species A evolve by 
synchronous updates in two sub-steps. Each site can be occupied by arbitrarily many 
particles. First the particles diffuse, i.e. all particles hop independently to a randomly 
chosen nearest neighbour. Then the pairs of particles react locally, either producing 
offspring or annihilating each other. To this end the population of n particles at a given 
site is divided into \n/2\ pairs, which branch independently with probability pi"-/ 2 ! anc [ 
annihilate otherwise. As a result, if n is very large, the fission process is exponentially 
suppressed, introducing effectively a soft constraint. 

Unlike restricted variants of the PCPD, the model introduced by Kockelkoren 
and Chate exhibits a surprisingly clean scaling behaviour, even when different order 
parameters are used (see figure 7). Since all reactions take place on single sites, the 
model seems to reach the scaling regime much faster than ordinary models with hard- 
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Figure 7. Soft-constraint model defined in [51] at criticality. From top to bottom 
are shown the decay of the particle density, the density of particles without solitary 
particles, the fraction of occupied sites and the fraction of occupied sites with at least 
two particles. A clean scaling behaviour is observed. After [51]. 

core exclusion. These results led the authors of [51] to the conclusion that the PCPD 
represents a single universality class different from DP and PC, characterized by the 
critical exponent 5 ~ 0.20. 

4-6. Surface critical behaviour 

While most of the attempts in understanding the critical behaviour of the PCPD 
concentrated on bulk quantities, some information of surface critical exponents is also 
available. It is a well-established fact, that in a semi-infinite system the critical behaviour 
near to the surface is in general different from the one deep in the bulk, see [83] for a 
recent review. For example, in the steady-state the order-parameter density should 
scale near criticality as pb ~ (p — Pc) 13 in the bulk but near the surface one expects 
Psurf ~ ijP — Pc) 131 where Pi is a surface critical exponent. The value of /?i may further 
depend on the boundary conditions. Here we shall consider 

(i) free boundary conditions, i.e. particles cannot cross the boundary. 

(ii) absorbing boundary conditions, i.e. particles at the boundary may leave the system 
with a rate D. In other words, the boundary reaction A —>■ is added. 

At criticality, one expects the finite-size scaling behaviour p sur f ~ L^ 131 /ly± and techniques 
to calculate the boundary density are available [37]. For the PCPD, Barkema and 
Carlon [54] have obtained estimates for the exponent ratio (3\/v± and we list their results 
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J 1 

model 


free 


absorbing 


method 


Rcl. 


PCPD 


0.72(1) 


1.11 


DMRG, D < 0.5 


[54] 


DP 


0.664(7) 




Monte Carlo 


[83] 




0.667(2) 




DMRG 


[37] 




0.6690(1) 




scries 


[84] 


PC 


0.73(1) 


1.11(1) 


Monte Carlo 


[83] 




0.720(2) 


1.10(1) 


DMRG 


[54] 



Table 5. Values of the surfrace exponent 0i/v± according to different methods for 
free and absorbing boundary conditions in some ID models. 
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Figure 8. PCPD interpreted as a cyclically coupled DP and annihilation process. 

in table 5 together with results for the DP and PC universality classes for comparison. 
In the PCPD, estimates only converge for D < 0.5. 

Taken at face value, the estimates of (3\/v±_ obtained for the PCPD appear actually 
to be compatible with those of the PC class but Barkema and Carlon carefully point 
out that their final estimates still depend on D, in particular for D large, and therefore 
may not yet give the correct asymptotic values. 

5. Related models 

5.1. Cyclically coupled spreading and annihilation 

In order to characterize the essential features of the PCPD from a different perspective, 
one of us (H.H.) suggested to interpret the PCPD as a cyclically coupled DP and 
annihilation process [40]. The idea is to separate the dynamics of pairs and solitary 
particles, associating them with two species of particles A and B. The A's, which 
stand for the pairs in the original PCPD, perform an ordinary DP process while the B's 
represent solitary particles that are subjected to an annihilating random walk. Both 
subsystems are cyclically coupled by particle transmutation, as sketched in figure 8. 

Such a cyclically coupled process was realized in [40] as a three-state model with 
random-sequential updates. Choosing particular rates it was shown that the process 
exhibits a phase transition with the same phenomenology as in the PCPD. Moreover, 
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Figure 9. (colour online) Typical spatio-temporal evolution of a binary spreading 
process starting from an initial seed. Particles of species A and B are represented as 
red and blue pixels, respectively. Plotting x/L 1 ^ 2 versus log 10 t the figure covers four 
decades in time. As can be seen, patches of high activity (red) are connected by lines 
of diffusing B-particles (blue) on all scales up to 10 6 time steps. 



the values of the effective critical exponents 

5 = 0.215(15), 2=1.75(5), ^ = 1.8(1) (37) 

were found to be in fair agreement with other estimates for the PCPD. Recent 
simulations for various values of D between 0.1 and 1 suggest that the exponents should 
be independent of the diffusion constant [85]. These results indicate that cyclically 
coupled processes following the react ion- diffusion scheme in figure 8 display the same 
type of transition as the PCPD. 

The separation into subprocesses illustrates an important feature of the PCPD, 
namely, the existence of two modes of spreading. As illustrated in figure 9, a 
typical spatio-temporal cluster is characterized by the interplay of a high-density mode 
dominated by self-reproducing and decaying A-particles, and a low-density mode of 
solitary diffusing 5-particles. A similar interplay of high- and low-density patches can 
be observed in the critical PCPD, see figure 1. 

Concerning the structure of a cluster a fundamental problem arises precisely at this 
point: As can be seen in figure 9, even after 10 6 time steps the solitary particles perform 
simple random walks over large distances. However, such a random walk is always 
characterized by the dynamic exponent z = 2, while all known simulations clearly 
indicate that the process as a whole spreads superdiffusively, i.e. z < 2. Therefore, 
the effective diffusion constant for solitary particles has to vary slightly under rescaling, 
meaning that a cluster such as in figure 9 cannot be scale-invariant. Therefore it seems 
that the process is still far away from the asymptotic scaling regime, even after 10 6 time 
steps. 
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r 


Pc 


z 


V 


5' 


5 = 0/v\\ 


v±. 


0.00 


0.04687(2) 


1.58(1) 


0.314(6) 


0.160(5) 


0.159(1) 


1.10(1) 


0.25 


0.05505(5) 


1.62(3) 


0.29(1) 


0.15(1) 


0.173(5) 


1.10(3) 


0.50 


0.06636(4) 


1.67(3) 


0.26(1) 


0.14(1) 


0.197(3) 


1.10(3) 


0.75 


0.08315(5) 


1.75(5) 


0.20(2) 


0.13(2) 


0.230(5) 


1.17(5) 


1.00 


0.1112(1) 


1.7(1) 


0.18(5) 


0.09(2) 


0.27(4) 


1.30(10) 



Table 6. Critical parameters of the generalized PCPD as defined in [47] for D = 0.1 
as a function of r. 

5.2. Interpolating between DP and PCPD 

A different explanation of the apparent non-universality in the PCPD was suggested 
by Noh and Park [47], who claim that the violations of scaling can be traced back to 
a long-term memory effect mediated by the diffusive background of single particles. To 
show this, they introduce a generalized version of the PCPD, as follows. Each site of the 
lattice can be either empty or occupied by a single particle and the following reactions 
are admitted 

AA -> 00 with rate (1 - D)p 

AA®, ®AA — > AAA with rate (1 - D)(l - p)/2 (38) 
A00 -> 0A0 with rate D 

A$A -> $AA with rate Dr (39) 

A$A -> 000 with rate D(l-r) 

which allows one to interpolate between the DP fixed point at r = and the critical 
behaviour of the PCPD at r = 1. In this model the branching probability depends 
on whether a collision is caused by a previous branching process or by diffusion. In 
addition, one restricts the calculation of averages to the sector where at least one pair 
(rather than two isolated particles) is present, in analogy with [49]. Since r = 1 in this 
model is not a special point, one may hope that corrections to scaling are more easy to 
control for < r < 1 and that the results could be extended to r = 1 at the end. 

Dynamic and static exponents are found from simulations at several values of r and 
some of the results of [47] are listed in table 6. As a consistency check, the exponent 
values are used to test the hyperscaling relation eq. (7) and a full agreement is observed. 

Comparison with table 3 shows that for r = 0, these results are in excellent 
agreement with the DP values while for r = 1, they fall within the range of values of the 
PCPD determined with different methods. In contrast to the limit D — * 0, where the 
exponents change abruptly, the exponents seem to vary smoothly when r is taken to zero. 
Based on these numerical observations Noh and Park [47] conclude that the transition in 
the PCPD may correspond to a line of fixed points with continuously varying exponents. 
As a possible explanation they suggest to consider the long-term memory imposed by 
solitary diffusing particles as a marginal perturbation of the underlying field-theory. 
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5.3. Parity- conserving PCPD 

As the PC class compared to DP exemplifies that an additional parity-conservation law 
may change the universality class of an absorbing phase transition the question arises 
whether the same happens in binary spreading processes. To answer this question in 
the context of the PCPD parity-conserving process 

2A -> AA , 2A -> . (40) 

was studied [42] . Surprisingly it turns out that parity-conservation does not change the 
nature of the phase transition, i.e., one observes the same phenomenological behaviour as 
in the PCPD with comparable effective critical exponents 5 ~ 0.23(1) and z = 1.80(5). 

To understand this observation, we note that there is another well-known example 
where parity conservation is irrelevant, namely, the annihilation process 2 A — > as 
compared to the coagulation process 2A — > A, which are known to belong to the same 
universality class. This is due to the fact that the even sector and the odd sector in 
the parity-conserving process 2A — > are essentially equivalent since in both of them 
the particle density decays algebraically until the system is trapped in an absorbing 
state (namely, the empty lattice or a state with a single diffusing particle). The present 
model is similar in so far as both sectors have an absorbing state. Contrarily, in the 
PC class only one of the sectors has an absorbing state, leading the authors of [42] to 
the conclusion that parity conservation changes universality whenever one of the two 
sectors has no absorbing state. 

5.4- Multicomponent binary spreading processes 

In order to search for further possible generalizations, Odor [46] investigated different 
variants of binary production-annihilation processes with several species of particles, in 
particular the reaction-diffusion processes 

(a) AB -> ABA, 2A -> 0, BA -> BAB, 2B -> 0. 

(b) 2A->2AB, 2,4 ^0, 2B -> 2BA, 25^0 

The process (a) seems to exhibit the same type of critical behaviour as the PCPD, wile 
in the process (b) the phase transition takes place at zero branching rate, where the 
exponents are different, e.g. f3 = 2. However, it would be misleading to conclude that 
the process (b) indeed represented a different universality class, rather the model is 
designed in a way that the nontrivial transition, which is possibly belonging to the same 
class as the PCPD, is shifted to an inaccessible region of the parameter space. Detecting 
a transition at zero branching rate does not necessarily mean that all realizations of the 
same reaction-diffusion scheme will yield such a trivial transition. 
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(a) Deposition (b) Evaporation (c) Hopping 

Figure 10. Dynamic rules of the PCPD-driven growth process: (a) Deposition of a 
dimer at the sites i and i + (b) Evaporation of a monomer from the right edge of a 
terrace at site i. (c) Hopping of a monomer between two adjacent islands, moving the 
hole between the islands in opposite direction. In all cases the spatially reflected rules 
(not shown here) have to be included as well. 



5.5. Roughening transition driven by the PCPD 

As demonstrated in [76] the PCPD may also play a role in a special class of models 
for interface growth which exhibit a roughening transition. The idea follows an earlier 
work by Alon et al. [86,87], where a deposition-evaporation model with a DP-related 
roughening transition was introduced. Later it was shown that by depositing dimers 
instead of monomers one obtains a different type of roughening transition which is 
related to the PC class [88,89]. 

As a key property of all these models, desorption is only allowed at the edges but 
not from the middle of deposited plateaus. Therefore the actual bottom layer of the 
interface, once it has been completed, cannot evaporate again which leads effectively to 
an absorbing phase transition at the bottom layer. Using this interpretation the pinned 
phase, where the interface is smooth, corresponds to the active phase of the underlying 
spreading process. However, if the growth rate is increased above a certain critical 
threshold the interface eventually detaches, meaning that the bottom layer enters an 
'absorbing' state. 

In order to study a PCPD-driven roughening transition one has to introduce 
appropriate dynamic rules that mimic a binary spreading process at the bottom layer. 
These rules involve three different physical processes, namely, deposition of dimers, 
surface diffusion, and evaporation of diffusing monomers (see figurelO). The dimers 
are deposited horizontally on pairs of sites at equal height, leading to the formation 
of islands. These islands are stable in the interior but unstable at the edges, where 
monomers are released at a certain rate. The released monomers diffuse on the surface 
until they either attach to another island or evaporate back into the gas phase. In the 
limit of a very high evaporation-rate a released monomer is most likely to evaporate 
unless it immediately attaches to an adjacent island at the next site, effectively moving 
the hole between the two islands in the opposite direction. 
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In order to explain how the growth model is related to the PCPD, let us consider the 
dynamic processes at the bottom layer (the spontaneously selected lowest height level of 
the interface). Interpreting sites at the bottom layer as particles A of a spreading process, 
the dynamic rules listed in figure 10 can be associated with certain reactions of the 
particles. For example, the deposition of a dimer corresponds to pairwise annihilation 
2A — > 0, while evaporation of a monomer can be viewed as the creation of a particle A. 
However, in the present model atoms can only evaporate at the edge of a terrace followed 
by two bottom layer sites, hence rule (b) in figure 10 has to be interpreted as a fission 
process 2A — > 3A. Otherwise, if there is only one such bottom layer site next to the 
edge, rule (c) applies which corresponds to a random walk of a single particle A. Thus 
the processes at the bottom layer resemble the dynamic rules of the PCPD. Clearly, this 
correspondence is not rigorous, especially because of second-layer nucleations. However, 
as in the case of DP- and PC-related growth models, the correspondence is expected to 
be valid asymptotically. 

In the case of DP-related growth processes, the behaviour at the first few layers 
above the bottom layer can be described by a unidirectionally coupled hierarchy of DP 
processes [90]. Similarly the dynamics at higher levels in the growth process in figure 10 
is expected to be related by a unidirectionally coupled sequence of PCPD's, which has 
not been investigated so far. 



6. Field-theoretic approaches 



Since the numerical simulations as a whole remain inconclusive, it is important to find 
and analyze a suitable field-theory that describes the phase transition in the PCPD. So 
far there are two important contributions in this direction. Some time ago Howard and 
Tauber [31] proposed a field-theory for the unrestricted PCPD, while recent attempts to 
devise a field-theory for the restricted case turned out to be unsuccessful [91]. Moreover, 
in a recent paper the pair correlation function in the unrestricted PCPD has been derived 
exactly by Paessens and Schiitz [53], leading to non-trivial and surprising results. The 
purpose of this section is to summarize the present state of knowledge. 

6.1. Field-theory of the unrestricted case 

Consider the unrestricted PCPD which is defined by the reaction-diffusion scheme 

2,4 ^0, lA^UA, 2A^3A, (41) 

where //, is, and a are the rates for annihilation, coagulation, and fission, respectively. 
In the continuum limit (see e.g. [92]) the classical master equation of this process 
corresponds to the field-theoretic action 



/d'x/d* 



<K<9t-W')0-Ml-0 )0 +^(1-0)CC (42) 
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where D is the diffusion rate while 0(x, t) and <p(x, t) can be thought of as particle 
creation and annihilation operators. Shifting the response field by = 1 + the action 
can also be written as 

S = Jd d x Jdt 0(<9 t -W 2 )0+(2/i+z/-a)00 2 +(/i+z/-2a)0 2 2 -a0 3 2 (43) 

Since higher-order terms will be generated under renormalization group (RG) 
transformations, it is convenient to write the action in the general form 



oo oo 

a d x / at \fi(dt - dv 2 )0 a p>o ( 44 ) 

J p=l q=2 

taking all possible reactions into account which are at least binary, i.e. quadratic in 0. 
In this notation the bare coupling constants of the unrestricted PCPD are given by 

Ai j2 = 2/j, + v - a , A 2i2 = li + v - 2<t , A 3i2 = —a (45) 

while all other bare coupling constants (but not necessarily the renormalized ones) 
vanish. Regarding the first two terms 00 2 and 2 2 a simple power-counting analysis 
yields the upper critical dimension 

4 = 2 (46) 

at which the dynamical critical exponent z — 2. In the annihilation phase the coupling 
constants and the fields carry the dimensions 

[A P ,,] = K 2 - {q ~ 1)d , [*] = **, = 1 (47) 

where k denotes an arbitrary momentum scale. Hence for given q all contributions with 
p running from 1 to oo are equally relevant and the usual renormalization group scheme 
in terms of a finite series of Feynman diagrams cannot be applied. At the critical point, 
however, the dimensions of the fields are given by [(f)] = [0] = n d l 2 [93] . 

In Ref. [31] the field-theory was found to be non-renormalizable so that no 
quantitative information on critical exponents can be obtained. Nevertheless it is 
possible to compute the critical point of the unrestricted PCPD exactly. To this end 
we note that the transition expressed in bare parameters takes place at 2fj, + v — a, 
where the cubic term 00 2 in the shifted action vanishes [31]. Since there is no way to 
generate this contribution under RG transformations to all orders, fluctuation effects do 
not influence the location of the transition, hence 

2[i + i/ = (T (48) 

determines the critical line to all orders of perturbation theory. In the inactive 
phase 2/j, + v < a annihilation and/or coagulation eventually dominate, leading to 
an asymptotic decay as t~ d l 2 for d < 2, while in the supercritical regime 2/x + v > a the 
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particle density grows rapidly and diverges within finite time. Consequently, there is no 
stationary state in the supercritical regime. 



6.2. Attempts towards a field-theory of the restricted PCPD 

Recently van Wijland, Tauber, and Deloubriere proposed and studied a field-theory for 
the restricted case. Their results, however, turned out to be unphysical [91], suggesting 
that a different field-theoretic approach is needed. 

Starting point of the theory by Wijland et al was to extend the field-theoretic 
action (42) by an additional mechanism which prevents the field from diverging in 
the supercritical regime. This mechanism may be implemented as a soft constraint by 
adding a reaction of the form 3A — > pA with p < 2. For example, the reaction 3A — ^ 
would lead to an additional contribution — r(l — 3 )0 3 in the action (42), corresponding 
to the bare coupling constants 

Ai, 3 = 3r, A 2 , 3 = 3r, A 3 , 3 = r . (49) 

Alternatively, the restriction may be incorporated directly in the field-theoretic action 
by inserting exponential damping factors, whose purpose is to 'switch off' the reactions 
when a certain particle density is exceeded [94]. To this end each reaction term is 
multiplied by an exponential function exp(— nv<j><f>), where n = max(p, q) is the maximal 
number of particles involved in the respective reaction and v is an additional coupling 
constant. When applied to the PCPD, one is led to the unshifted action 

S = J d d x J dt 4>{d t - DV 2 )0 - fi(l - 2 )0 2 e" 2 ^ (50) 
- v(l - 0U0 2 e~ 2 ^ + a(l - 6)6 2 6 2 e- 3v ^ 



Roughly speaking, the constant v determines the field amplitude needed to declare 
location in space-time (representing a site) as occupied. Since the argument of the 
exponential function has to be dimensionless and [00] = n d the coupling constant v is 
irrelevant, allowing the exponential functions to be expanded. To first order in v this 
gives rise to the additional non-zero bare coupling constants 

Ai,3 = v(3a -Aii- 2v) 
A-2,3 = v(9a — 6/i — Av) 
A 3 ,3 = v{9a -2fi- 2v) 
A 4i3 = 3av 

while higher-order contributions lead to additional non-vanishing bare coupling 
constants of the form A p>q oc v q ~ 2 . 

In the inactive phase all terms of the leading series A Pi 2 are marginal and thus have 
to be retained while all terms of higher order are irrelevant. In order to handle the series 
of marginal contributions, Wijland et al studied the corresponding generating function 
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G(x) = Yl^Li ^p,2(p 2 4> p , recasting the flow equations of all coupling constants A p 2 into a 
single functional renormalization group equation for the renormalized counterpart of G. 
Alternatively, one may work at the critical point, where the naive field dimensions of <j) 
and <f> are equal, leading to a renormalizable field-theory [93]. Although a fixed point was 
shown to exist, it turned out to be unreachable from any physically meaningful initial 
condition. In addition, the fixed point itself is characterized by a negative variance of 
the order parameter fluctuations, which would have been unphysical as well. Finally, the 
critical exponents, in particular the exact result z = 2, would have been incompatible 
with the existing numerical estimates which clearly indicate that z < 2. 

The failure suggest that a field-theoretic action in terms of a single coarse grained 
field (and its associated response field 0) may be incapable to describe the properties 
of the restricted PCPD. In fact, a single field for the density of particles may fail to 
describe the complex spatio-temporal interplay of particles and pairs, questioning the 
validity of the continuum limit on which Eq. (50) is based. Instead it may be more 
promising to devise a field-theory based on two different fields for single particles and 
pairs, separating the two dynamical modes in a similar way as in the cyclically coupled 
models of Sec. 5.1. However, a consistent formulation of such a field-theory is not yet 
known. 

6. 3. On the impossibility of a Langevin equation for the PCPD 

Langevin equations are very popular in non-equilibrium statistical physics because of 
their intuitive simplicity and their enormous success in well-known systems such as 
DP. Therefore it has become customary to deal with phenomenological or even guessed 
Langevin equations even when a rigorous derivation is not available. For example, in 
the context of the PCPD several authors [31,51,55,95] have been discussing a Langevin 
equation of the form 

^p(x, t) = (a- 2p)p 2 (x, t) - cp 3 (x, t) + DV 2 p(x, t) + p(x, t)£(x, t) , (51) 

where £(x,t) is a white Gaussian noise with zero mean and the correlations 

(t(x, t)£&, 0> = (4a - 2/i) (x - x')5(t - . (52) 

This Langevin equation, if indeed realized, would be particularly appealing since it can 
be related to the problem of non-equilibrium wetting [55]. However, as a mesoscopic 
description a Langevin equation assumes that all features of the process, in particular 
all particle-particle correlations, could be captured in terms of a single coarse-grained 
particle density. This works only in few cases (e.g. in DP), but not in the case of the 
PCPD, where anticorrelations produced by the annihilation/coagulation process play 
an important role. 
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In this context it is useful to recall that a Langevin equation can be derived from 
the field-theoretic action by transforming the terms of the form 2 9 into a noise and 
computing the functional derivative with respect to the response field 0. This procedure 
requires that terms involving cubic or higher-order powers of the response field are fully 
irrelevant. However, in the case of the PCPD cubic terms of the form 3 9 have to be 
retained in the action. Therefore, any attempt to describe the PCPD by a Langevin 
equation misses the physics of these terms and thus should lead to wrong results. In 
addition, one would have to add an additional term £d{&, t) on the r.h.s. of the Langevin 
equation accounting for diffusive noise with the correlations 

{Sd(x, t)£ D (x', 0) = -D p(x, t) V 2 5^(x - x')S(t - t') (53) 

which is expected to play an important role in the PCPD. However, it is at present not 
clear what consequences the inclusion such an additional term would have [93]. 

6.4- Pair correlation functions in the unrestricted PCPD at criticality 

A different approach has been followed by Paessens and Schiitz [53] and leads to a non- 
trivial, unexpected result. They consider the unrestricted version of a generalized PCPD 
with the reactions 

mA — > (m + k)A with rate v 

pA — > (p — l)A with rate \i (54) 

<-> with rate D . 

Models of this kind may be thought of as realizing a bosonic field-theory such as in 
Ref. [31] and are specified by four integers (m,p,k,l) and three rates u, fi and D. 
Paessens and Schiitz start from the master equation (19) and rewrite the quantum 
Hamiltonian H in terms of space-time-dependent creation and annihilation operators 
a\x,t) and a(x,t) such that the particle operator n(x,t) = (x,t)a(x,i). In what 
follows, the averages 

(n(x,t)) = (a(x,t)) , (n(x,t) 2 ) = (a(x,t) 2 ) + (a(x,t)) , (55) 

will be considered. Specializing to the case m = p, the critical line between the 
absorbing phase and the active phase (with an infinite particle density) is located at 
\i c = vk/l. For /i = /i c , the spatially averaged particle density J^da; (a(x,t)) = p 
is a constant, where Q is a spatial domain with volume Remarkably, it turns out 
that at criticality and for m = p, the system of equations of motion for the variables 
{(a(x,t)), (a(x,t)a(y,t))} closes if either m = 1 or m = 2 [53]. 

The case m = 1 is a 'bosonic' version of the ordinary contact process and has 
been considered as a model for biological clustering [63,96]. It can be shown that the 
fluctuations (a(x,t) 2 ) of the particle density diverge for t — > oo in dimensions d < 2, 
while they remain finite for d > 2 [53,63]. 
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The PCPD is described by the case m = 2. After a time rescaling t i— > t/(2D), and 
assuming translation-invariant initial conditions, consider the pair-correlator 

F(r, t) := (a(x, t)a(x + r, t)) (56) 

(which is independent of x) and the reduced coupling 

v k(k + I) , r _, 

" = £>~2 (57) 



Then, under the stated conditions, the following equation holds on a hypercubic lattice 

d d 

-F(r, t) = J2 i F (r ~ e h t) + F(r + e h t) - 2F(r, t)} + 5 r , aF(0, t) (58) 
i=i 

where is the unit vector in the i th direction. If as initial condition one takes a Poisson 
distribution F(r, 0) = pi, one obtains the following Volterra integral equation 

F(r, t)=p 2 + a [ dr b(r, t - r)F(0, r) (59) 
Jo 

b(r,t) =e- 2dt I ri (2t)---I r . d (2t) (60) 

where I r is a modified Bessel function of order r. For r = 0, the same equation describes 
the kinetics of the mean spherical model, see [97-99] but it is not yet clear if there is a 
deeper reason for this relationship or whether it merely occurs by accident. The exact 
solution leads to the following long-time behaviour, at p = p c [53] 

{exp(t/r) ; if a > a c 
t v ; iia = a c (61) 

Foo ; if oi < a c 

where a c is a known critical value such that a c = for d < 2 but it is finite for d > 2, 
in addition r = Cd ( a ~ ac ) the exponent v = d/2 — 1 for 2 < d < 4 and v — 1 for 
d > 4, = Po(l ~ 0L l 0L c) X an d ^ is a known constant. As can be seen the two-point 
correlator diverges for large times if the reduced rate a is large enough. This result 
can be generalized to yield F(r, i) exactly, with a similar conclusion. In particular, one 
obtains the dynamical exponent z = 2 [53]. 

Since a measures the relative importance between the reaction and the diffusion 
rates, the result (61) means that for d > 2 at high values of D the distribution of 
the particle density along the critical line is relatively smooth and changes through 
a 'clustering transition' (which occurs at a tricritical point) to a rough distribution 
for small values of D. This transition manifests itself in the variance c(t) 2 : = 
(n(x,t) 2 ) — (n(x,t)) 2 and not in the mean density {n(x,t)) which is purely diffusive 
and has a constant average. A convenient order parameter of the clustering transition 
is F^ 1 [53]. On the other hand, for d < 2 there is a single transition along the critical 
line. 
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These exact results for the bosonic PCPD [53] surprisingly agree with the 
conclusions of pair mean-field theory and with the numerical results of [39] which 
suggested the presence of two distinct transitions along the critical line. That scenario 
had, in the light of more recent simulations in ID, somewhat fallen into disfavour, see 
section 3. In fact, in hindsight one observes that the existing cluster approximations of 
the PCPD all explicitly restrict to clusters of a linear shape. At the level of the pair 
approximation, at most two-site clusters need be considered and there is no qualitative 
difference between one and several space dimensions. On the other hand, beginning with 
the triplet approximation, the geometrical shape of the clusters also becomes important 
and one will have to distinguish between linear clusters such as «o« or •••o (in ID only 

such clusters occur) and non-linear ones as or which are meant to indicate 

o« • o 

more complicated geometrical forms as they will only occur for d > 2. Indeed, truly 
three-dimensional clusters arise for the first time in the quartet approximation (N — 4). 
The triplet and higher approximations of the form as discussed in section 3 only take 
linear clusters into account and should hence be expected to be more adapted to a truly 
one-dimensional system than the pair approximation. Then a single transition along the 
critical line should be expected for the N > 3 approximation while the different result 
of the pair approximation might be related to the existence of two transitions for d > 2. 

It would be very interesting to see whether two transitions exist along the critical 
line for the restricted PCPD in d > 3 dimensions. The interaction terms in existing 
field-theory studies might be too close to the site-approximation to be able to see this 
and a study based on the master equation is called for. 

7. Possible generalizations 

7.1. Higher-order processes 

As shown in the preceding sections, the unusual critical behaviour of the PCPD is 
related to the fact that only pairs of particles can react while individual particle form 
a diffusing background. Therefore, as an obvious generalization, it is near at hand to 
study higher-order processes such as 

nA -> (n + I) A , nA -> . (n > 3) (62) 

Without diffusion, the number of absorbing steady-states increases exponentially with 
the lattice size L. For example, for the triplet process (n = 3) and in one dimension, we 
find 

N s UL)cN o9 L , N =[ 1 ; P e ™ dicb - C - (63 ) 

staU ; oy ' 1 1.137... ; freeb.c. V ; 
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Ref. 


5 


z 


(3 


V\\ 


Park et al. [95] 


0.32(1) 


1.75(10) 




2.5(2) 


Kockelkoren and Chate [51] 


0.27(1) 


1.8(1) 


0.90(5) 




Odor [45] 


0.33(1) 




0.95(5) 





Table 7. Estimates for the critical exponents of the ID triplet process TCPD. 



absorbing steady-states, where g ~ 1.839 . . . (see appendix A). On the other hand, with 
single-particle diffusion, only n absorbing steady-states remain. 

The central question is whether such a system still exhibits a non-trivial phase 
transition. As a first naive approach, the simple mean-field equation 

^p(x, t) = bp n (x, t) - cp n+1 (x, t) + DV 2 p(x, t) (64) 
predicts the mean-field exponents 

oMF i ,, MF / , , MF „ / C r\ 

p — 1, v± = n/2, z/|) = n. (65) 

in sufficiently high space dimensions, including DP (n = 1) and the PCPD {n = 2) 
as special cases. The question posed here is whether the higher-order processes, in 
particular the diffusive triplet-contact process (TCPD) with n — 3, still displays a non- 
trivial critical behaviour in one spatial dimension. 

The TCPD was introduced in [95], where numerical simulations seemed to indicate 
a possible non-trivial behaviour at the transition. This point of view was confirmed 
by Kockelkoren and Chate [51], who established the ID TCPD as an independent 
universality class and with an upper critical dimension d c — 1. However, Odor [100] 
questioned this claim by showing that the simulation results are compatible with mean- 
field exponents combined with appropriate logarithmic corrections (see table 7). 

Even more confusing are the findings for the quadruplet contact process n — 4. 
In [51] the expected mean- field exponents are obtained while Odor finds numerical 
evidence for the presence of fluctuation effects [100] . The origin of these discrepancies is 
not yet understood. However, when comparing the results one should keep in mind that 
the model used by Kockelkoren and Chate is particularly suitable to simulate higher- 
order processes. In their model the reactions are local for any n, while in models with 
hard-core exclusion extended strings of n particles have to line up before they can react. 

7.2. Towards a general classification scheme 

Extending these studies Kockelkoren and Chate [51] suggested that transitions in a 
react ion- diffusion processes of the form 



mA -> (to + k)A , nA -> (n - I) A , 



(66) 
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m\n 


1 


2 


3 


4 


1 


DP 


DP/PC 


DP 


DP 


2 


DP 


PCPD 


PCPD 


PCPD 


3 


DP 


DP 


TCPD 


TCPD 


4 


DP 


DP 


DP 


? 



Table 8. General classification scheme of absorbing phase transitions in the ID models 
(66) as proposed by Kockelkoren and Chate [51] and modified after [102]. 

can be categorized according the the orders m, n of creation and removal, while the 
numbers k and I determine additional symmetries such as parity conservation. In terms 
of m, n, they propose the general classification scheme shown in table 8. In addition they 
investigate the role of parity-conservation in more detail, finding that such a symmetry 
does not alter the universality class whenever every sector includes an absorbing state. 
This conjecture is in agreement with the observation of Ref. [42] that an additional 
parity-conservation does not change the critical behaviour of the PCPD (see [101] for 
recent evidence of a PCPD transition in the diffusive 2A — > 3A, 4A — > system, where 
777, — 2,71 = 4). 

8. Summary and Outlook 

In spite of intensive and prolonged efforts, the critical behaviour of the one-dimensional 
PCPD is not yet understood. It is surprising to see that the field-theoretical and 
numerical methods used to study the PCPD-transition have led to such wildly different 
conclusions as listed in section 1. It appears that almost the only piece of information for 
which different methods obtain the same result is the location p c {D) of the critical point 
(if applicable at all). This continuing discrepancy is all the more astonishing because 
the same techniques yield nicely consistent results when applied to models in the DP or 
PC universality classes. By itself, this observation is a clear indication of a subtlety in 
the behaviour of the PCPD which is not present in those other models. 

Specifically, at present there seems to be some majority opinion, based on numerical 
work, in favour of a single new universality class along the critical line and characterized 
by critical exponent values of the order 

5^0.2 , zwl.7-1.8<2 , r/^0.5 , (3/v ± « 0.5 (67) 

On the other hand, a recent analytical study [53] of the bosonic field-theory finds a 
dynamical exponent z = 2 and a single universality class along the critical line, at least 
in one dimension. Their value of z = 2 clearly is in disagreement with the numerical 
results from fermionic models. At the time of writing, it is not clear how this should be 
explained. 
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Besides the majority opinion, several alternative scenarios have been proposed. 
These include the possibility of more than one transition along the critical line (which 
is known to be exact for the unrestricted bosonic version in d > 2 dimensions), 
continuously varying exponents and even the cross-over to directed percolation after 
very long times. The last possibility would mean that the diffusive background of 
solitary particles would become irrelevant in the limit t — > oo which goes against the 
current understanding of the bosonic field-theory. 

All in all, we are not (yet) able to decide between the several scenarios proposed. 
Pending further insight, the reader might find some comfort in a quotation of 'Thales' 

Sei ruhig - es war nur gedacht. J. W. Goethe, Faust II (1830) 

or alternatively with 

It were not best that we should all think alike; it is the difference of opinion 
which makes horseraces. M. Twain, Pudd'nhead Wilson (1894) 

We look forward to future insightful studies shedding more light on the so simple-looking 
PCPD. 
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Appendix A. The absorbing states in the pair-contact process 

Consider the ID pair-contact process on a lattice with L sites. We derive the number 
of absorbing steady-states iV stat (L) given in eq. (1), following [19]. The main ingredient 
is that any state which only contains isolated particles is an absorbing steady-state. 

First, we consider free boundary conditions. Obviously, iV stat (l) = 2 and iV stat (2) = 
3. Now fix the leftmost site. If it is occupied, its neighbour must be empty in order 
to obtain a steady-state and it remains to consider an open chain of L — 2 sites. On 
the other hand, if the leftmost site is empty, one considers the remaining open chain of 
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L — 1 sites. We thus have the recursion 



N stat (L) = A stat (L - 2) + A stat (L - 1) 



(Al) 



Because of the initial conditions, 



A stat (L) =F L+1 




(A2) 



9+~ 9- 



is the (L + l) th Fibonacci number, where g± — (1 ± y / 5)/2. 

On the other hand, for periodic boundary conditions, we fix one of the sites. If that 
site is occupied, both its left and right nearest neighbours must be empty and we are 
left with an open chain of L — 3 sites. But if that site is empty, we are left with an open 
chain of L — 1 sites. Therefore 



From these results we recover (1) for L sufficiently large. 

In general, N stat (L) is found from its generating function Q(z) = J2T=o-^stat{L)z L . 
For the kind of application at hand, Q(z) = z(2 + z)/(l — z — z 2 ) is a simple rational 
function of z and its coefficients can be obtained in a closed form as follows. 

Lemma: Let P(z) and Q(z) be polynomials of orders p and q, respectively and letp < q. 
If one has in addition Q(z) = qo(l — zpi) ■ ■ ■ (l — zp q ) such that the pi are pairwise distint 
for i = 1, . . . , q and q is a constant, then 



In [103] this is proven in an elementary way, here we give a function-theory proof 
which easily generalizes to the case when Q(z) has multiple zeroes. We have 



where we have set u = 1/w. The contour C is a circle around the origin with a radius 
smaller than the convergence radius of f(z) and C encloses all the points pj, j — 1, . . . , q. 
Since p < q, there is no singularity at u = and the only singularities of u^ 1 f{l/u) are 
the simple poles located at u — pj. If we concentrate on the singularity at u — pi, we 
have uQ(l/u) = (u — pi)qo YYj=2 (1 ~ Pj/ U ) suc h that the product is regular at u — p\. 
On the other hand, the derivative Q'(l/p\) = — piqoY[ Q j=2 (1 ~~ Pj/Pi)- Therefore, close 
to u ~ pi, we have 



Aff (£) = A stat (L - 3) + A stat (L + £ 



(A3) 




(A4) 




uQ(l/u) 



(u - Pi) 



Q'(l/Pi) 
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From the residue theorem and summing over all simple poles of f(z) the assertion 
follows. q.e.d. 

Other processes can be treated similarly. As an example, consider the triplet- 
contact process without diffusion as introduced in section 7.1. Absorbing states are 
those without any triplet of neighbouring occupied sites. On a ID lattice with free 
boundary conditions and L sites, we have the recursion 

iV stat (L) = N stat (L - 1) + iV stat (L-2) + iV s tat (L - 3) (A5) 

together with the initial conditions As tat (l) — 2, N stat (2) = 4 and iV B t a t(3) = 7. 
Therefore the generating function is Q{z) = z(2 + 2z + z 2 )/(l — z — z 2 — z 3 ) and its 
coefficients are easily read off. For periodic boundary conditions, we find 

N^\L) = N stat (L - 1) + N stat (L - 3) + 2A stat (L - 4) (A6) 

For L sufficiently large, we thus arrive at eq. (63). 



Appendix B. On finite-size scaling techniques 

In systems like the PCPD, where the the critical line separates a non-critical ordered 
phase from a critical disordered phase, the analysis of finite-size data needs some 
modifications with respect to the usual situation, where both the ordered and the 
disordered phases are non-critical. To be specific, we shall discuss here the finite- 
size scaling of the lowest gap T of the quantum Hamiltonian H, see section 4, and 
following [19]. As usual, we begin with the finite- size scaling form 

r L ( P ) = L-f((p-p e )L 1 ^) (Bl) 

where p is the control variable which measures the distance to the critical point p c 
(the dependence on D is suppressed throughout) and / is assumed to be continously 
differentiable. One expects the asymptotic behaviour, see eq. (23) 



where a = a{p) is a constant. In the usual case, one would have instead T L {p) ~ Y 
for p > p c . From (B2), one finds for the scaling function 



j exp(-A|3|^) ; if 3 -> -oo 

;if 3 -+oo (B3) 



where A is a positive constant. Therefore, since f(i) > 0, it follows that for z < 2 the 
scaling function / must have a maximum at some finite value 3 max - Next, the logarithmic 




Figure 11. Plot of the function Yl(p) in the ID PCPD as a function of p for several 
lattice size L and D = 0.5. The inset shows the case D = 0.2. After [19]. 



derivative (24) becomes 

Y _ : . ln[/( 3+ )//(j-)l 

where 3± = (p — p c )(L ± l) 1 /"- 1 -. Furthermore, it is easy to see that in the scaling limit 
p — > p c and L — > oo such that 3 = (p — p c )L l ^ u ^ is kept fixed, one has 

lim dp -\ LV- (z _ 2)3 -i • if^ +00 

Provided that z < 2 and < 2, there must exist a finite 3* such that dli/dpl » = 0. 
However, since 

that maximum value of Yl cannot be used to estimate the dynamical exponent z [19]. 
Rather, one has to form first a sequence of estimates pi of the critical point p c from 
the above extremum criterion which should converge according to pl — p c + 3*L _1 ^ X . 
Having found p c , estimates of z can finally be obtained from eq. (25). 

We finally point out that the habitual method of intersecting two curves Y^{p) and 
Yl'(p) for two lattice sizes L, L' does not work in general. For example, it is known for 
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the ID PCPD with free boundary conditions, that the curves Y L {p) do not intersect for 
different values of L [19]. We illustrate this in figure 11. 



Appendix C. On universal amplitudes in nonequilibrium criticality 

We recall the phenomenological scaling arguments about universal scaling amplitudes 
in non-equilibrium critical points and derive eqs. (33,34) for sytems in the directed 
percolation universality class, following [44]. We shall denote time by r\\ and space by 
r±. The distance from the steady-state critical point is mesured by t and h denotes an 
external field (e.g., for directed percolation t = p — p c and h is the rate of a process 
— > A). We shall assume translation invariance throughout. The arguments involves 
tracing the non-universal metric factors, generalizing a similar by now classic line of 
thought for equilibirium systems [80]. This is most conveniently done in the bulk first 
before finally introducing lattices of finite extent. 

Physical quantities of interest are the mean particle density p, the survival 
probability P and the pair connectedness function G = G(r' ± , r'p r±, r\\), which is defined 
as the probability that the sites (r' ± , rf,) and (r±, ry) are connected by a direct path [1,2]. 
Because of translation invariance G = G(r' ± — r±, r|| — ry). These quantities are expected 
to satisfy the scaling behaviour 

p(r ± , rf t, h) = b-*p (^, W, hff^) = D lp ^£ ± (g, D ^; D 2 h\t\~^ 

P(r ± ,r f t,h) = b-*-P (^-, r -±-W\hb^) = D lP ^r P ^ (j^,D ^,D 2 h\t\-y^ 

G(r ± ,r f t,h) = b~*°G = D ^; D 2 h\t\-^ 

(CI) 

where the x's are scaling dimensions and y t> h renormalization group eigenvalues (in 
particular y t = l/v±), the D's are non-universal metric factors, £,T,Q are universal 
scaling functions where the index distinguishes between the cases t > and t < 0, 
= ^o|^| _iy± is the spatial, £y = £±/ D is the temporal correlation length and z is the 
dynamical exponent. 

In the steady state, and for h = 0, one expects p ~ t 13 and P ~ t^' . In general, the 
two exponents j3 and f3' are distinct from each other. For spatial translation invariance, 
the dependence on r± drops out for both p and P and in the steady state (i.e. r\\ — > oo) 
one has 

p(t,h) = D lp £-^£ ± (D 2 h\t\~y^) \t\f> 

P(t,h) = D XP i~ 0lv ^T± {D 2 h\t\- yh ' yt ) \tf (C2) 
where x p = f3/v±, xp = /3'/v± and £ ± = \im ri ^ OQ £ ± and similarly for T. We also 
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consider the auto-connectedness (that is r± = r' ± ) in the steady state 

G{0,oo;t,h) =: G(t,h) = D 1P £ XG g ± (D 2 h\t\- yh/yt ) \t\ XGU ^ (C3) 

In the active phase (t > 0), the surviving clusters will create an average density ~ \t\^ 
in the interior of the spreading cone. Therefore, the auto-connectedness should in the 
steady state saturate at the value [3] 

G(t,h) = p(t,h)P(t,h) (C4) 

Comparison of the scaling forms then yields, setting h — 0, 

x G = ((3 + P')/u ± , Dig = D lp D 1P ^ ^ ) (C5) 

Usually, xq = d — 6z is expressed in terms of the initial critical slip exponent 9 [104], 
which makes it apparent that the expression (C5) is in fact a generalized hyperscaling 
relation, see [3]. 

Next, we consider the total mass M of the cluster, given by 

M(t,h):= [ dV rdr ]] G(r ±1 r r ,t 1 h) = ^^ v± ^(D 2 h\t\^) (C6) 
Jn d Jo u o 

where eq. (CI) was used and is a new universal function related to Q ± . Also 

1 = du ± + u\\ -P-(? (C7) 

which is the analogue of the hyperscaling relation of equilibrium systems. 

While the discussion so far has been completely general, we now appeal to two 
properties which are valid for systems in the directed percolation universality class, but 
need not be generically valid. First, we consider a directed percolation process in the 
presence of a weak field h (physically, h parametrises the rate of a particle creation 
process — > A). A site at a given time becomes active if it is connected with at least 
one active site in the past, where a particle was created by the field. The number of 
such sites is equal to the cluster size, the probability to become active is given by the 
density 

p(t, h) ~ 1 - (1 - h) M(t ' h) ~ hM(t, h) (C8) 
for h small. Therefore, 

dp(t,h) 



M(t,0) 



(C9) 

h=0 



dh 

Comparison with the scaling forms for p and M leads to 

Vh/Vt = (3 + 7 , D 1P = D D 2 ^ +l)/l/± A ± (CIO) 
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where A ± is an universal amplitude. Second, directed percolation is special in the sense 
that there is a 'duality' symmetry which can be used to show that [23] 

p(t,h) = P(t,h) (Cll) 

As a consequence, f3 — (5' and D lp = D 1P for directed percolation and we thus have, 
combining eqs. (C2,C6,C10) 

p(t, h) = D D 2 £ d - z \tfM± (D 2 h\t\~^) 

M(t, h) = D D 2 2 ^ d ~ z \t\-' f Mi (D 2 h\t\-^) (C12) 

with universal functions Ai^(x) = d n Ai ± (x)/dx n and where the hyperscaling relation 
eq. (C7) has been used. Finally, we define a new function p = pit, h) by p(t, h) = 
dp(t,h)/dh, which implies 

pit, h) = D ^ d - Z \t\^ d+Z) ^M ± (D 2 h\t\-^) (C13) 

In particular we have because of £y = £±/D that 

(i(t,0)Zi(t, 0)f|| (*, 0) -> univ. constant, (C14) 

again quite analogous to an equilibrium result, see [80]. 

At last, we consider a geometry of finite size L in space but of infinite extent in 
time. As usual in finite-size scaling [80], we postulate that in this finite geometry merely 
the scaling functions are modified 

M$ = M$(D 2 h\t\-^I4Z x ) (C15) 

and without introducing any further metric factor. Indeed, we can then scale out L and, 
because of eq. (C14), arrive at the scaling forms (33). The universality of the moment 
ratio (34) is obtained from the analogous extension of (CI) as follows 

L {L) ^±Io^nrlG(r ± , n ;L/^) 



(rl) 



L(L) ddr -L Jo°° dr ii g (?-l> i ; L l £0 

, n L(L/t ± ) ddr± io°° dr ll r T XG G ± ( r ±, r\\;L/Z A 



W) d ^ /o°° dr ll rToQHrj., r\\;L/U) 

= L n Z n (L/t ± ) (C16) 

where Q(L) is a ci-dimensional hypercube of linear extent L and S n and S n are universal 
functions. 
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